{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Weighted Linear Fit\n", "\n", "*June 9, 2021*\n", "\n", "In this notebook we will fit a linear function to a set of experimental data. The fit will be weighted by the measurement uncertainties.\n", "\n", "Updated by Jordan Andrews on June 9, 2021 with the use of np.polynomial.Polynomial([a, b, c]).\n", "\n", "In this example, our data will be the voltage across and the current through a resistor. The slope of our data will tell us the resistance of the resistor. First import the *NumPy* module." ] }, { "cell_type": "code", "execution_count": 72, "metadata": {}, "outputs": [], "source": [ "import numpy as np" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Enter the data as arrays." ] }, { "cell_type": "code", "execution_count": 73, "metadata": {}, "outputs": [], "source": [ "Irms= np.array([0.00907, 0.00794, 0.00642, 0.00472, 0.00296, 0.001911, 0.000467, 0.0000469])\n", "errIrms= np.array([1.9e-4, 1.8e-4, 1.6e-4, 1.5e-4, 1.3e-4, 1e-4, 5e-5, 1e-5])\n", "Vp2p = np.array([5.62, 4.94, 4.02, 2.96, 1.86, 1.19, 0.302, 0.047])\n", "errVp2p = np.array([0.04, 0.04, 0.04, 0.04, 0.04, 0.03, 0.01, 0.008])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "When the data was collected the current was measured as rms values, convert the peak-to-peak voltage measurements to rms by dividing by $2\\sqrt{2}$." ] }, { "cell_type": "code", "execution_count": 74, "metadata": {}, "outputs": [], "source": [ "Vrms = Vp2p/(2*np.sqrt(2))\n", "errVrms = errVp2p/(2*np.sqrt(2))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Notice that the current data has larger relative errors than the voltage data (typically by a factor of two or more):" ] }, { "cell_type": "code", "execution_count": 75, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([2.9432194 , 2.79974811, 2.5046729 , 2.35169492, 2.04222973,\n", " 2.07570208, 3.23340471, 1.25266525])" ] }, "execution_count": 75, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Irel = errIrms/Irms\n", "Vrel = errVrms/Vrms\n", "Irel/Vrel" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ " Because the relative error in current is larger, we will plot the data as current vs voltage and show the error bars only along the y-axis. Use plt.errorbar(x,y,e) from the *Matplotlib* module." ] }, { "cell_type": "code", "execution_count": 76, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEGCAYAAABy53LJAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAAfVUlEQVR4nO3dfZwcVZ3v8c+XmeGigQBKojiDBtcAQsAEZxlQ2QeVhGETwogPoBd8oSsX7rK7IGsEVkgLu4igV8UHuFk3K3A1iBeCCTJOgFXU1QAJGUgiJGbxgQm4xL0sDwFlZvjdP6qGNJ2Znqpkaroz/X2/Xv2arqpzqn/dr4JfTp1T5ygiMDMzy2q3WgdgZma7FicOMzPLxYnDzMxyceIwM7NcnDjMzCyX5loHMB7222+/mDZtWq3DMDPbpaxevfp3ETGlcn9DJI5p06axatWqWodhZrZLkfTr4fb7VpWZmeXixGFmZrk4cZiZWS5OHGZmlosTh5mZ5eLEYWZmuThxmJlZLk4cZmaWixOHmdkEUSqVkDTqq1Qq7dTnqBEWcmpvbw8/OW5mjWRwEKZMOZ3+/hksWbKAzk5oasp3DkmrI6K9cr9bHGZmE8zgIHR1bWXq1HM5++xBFi7cQFfXVgYHx+b8DTFXlZlZI+nuhs2b+1i7toOWlgH6+y+ho2M93d0HMXfuzp/fLQ4zswlmzRqYPftWWloGAGhpGWDOnKX09o7N+Z04zMwmmFmzYMWKk+jvT24q9fc309PTxcyZY3N+Jw4zswmmsxNaW9uYMeNeFiy4nI6O9bS1tdLZOTbndx+HmdkE09QES5dOYsqUL3DNNYexZMlBOzSqaiRucZiZTRDlz3E0N4snn7yBZ5+9gHnzkm0/x5GDn+MwM8vPz3GYmdmYcOIwM7NcnDjMzCwXJw4zM8vFicPMzHJx4jAzs1ycOMzMLBcnDjMzy8WJw8zMcnHiMDOzXJw4zMwsFycOMzPLxYnDzMxyceIwM7NcCk0cko6XtEHSJkkXDHNckq5Ojz8o6cjR6kqaKWmlpF5JqyQdVeR3MDOzlysscUhqAr4KdAKHAqdKOrSiWCcwPX2dCVyToe6VwKcjYiZwSbptZjZuyhdMqvba2QWT6lWRLY6jgE0R8UhEvADcCMyvKDMfuD4SK4F9JO0/St0AJqfv9wYeK/A7mJltp1QqERFEBAMDwb77nsaee36W5cuT7aFjThz5tQKPlm33pfuylKlW91zgKkmPAp8DLhzuwyWdmd7KWrVly5Yd/Q5mZiMaHISurq1MnXouZ589yMKFG+jq2srgYK0jK1aRiUPD7Ktcp3akMtXqng2cFxEHAOcB/zzch0fEoohoj4j2KVOmZAzZzCy77m7YvLmPtWs7uPLKi1i5cgZ9fZvp7q51ZMUqMnH0AQeUbbex/W2lkcpUq/th4Jb0/XdIbmuZmY27NWtg9uxbaWkZAKClZYA5c5bS21vbuIpWZOK4D5gu6UBJuwOnAMsqyiwDTk9HVx0NPBURj49S9zHgT9P37wR+UeB3MDMb0axZsGLFSfT3NwPQ399MT08XM2fWNq6iNRd14ogYkHQO0AM0AYsjYr2ks9Lj1wK3AycAm4DngDOq1U1P/THgS5Kagd+TjMYyMxt3nZ2waFEbM2bcy/z53+fOO0+mra2Vzs5aR1YsRVR2O0w87e3tsWrVqlqHYWYTRKlU4tOf/nS6tRvJkwMzgV6gG3gRgIULF+7SI6skrY6I9u32O3GYmdlwRkocnnLEzMxyceIwM7NcnDjMzCwXJw4zM8vFicPMzHJx4jAzs1ycOMzMLBcnDjMzy8WJw8zMcnHiMDOzXJw4zMwsFycOMzPLxYnDzMxyceIwM7NcnDjMzCwXJw4zM8vFicPMzHKpuua4pD2AucCxwOuA54F1wPfK1gA3M7MGMmLikFQC5gE/BO4BngD2AA4CrkiTyvkR8WDxYZqZWb2o1uK4LyJKIxz7X5KmAq8f+5DMzKyejZg4IuJ7lfsk7QbsGRFPR8QTJK0QMzNrIKN2jkv6lqTJkiYBPwc2SPpE8aGZmVk9yjKq6tCIeBo4Cbid5PbUaUUGZWZm9StL4miR1EKSOL4bEf1AFBqVmZnVrSyJ438DvwImAT+S9Abg6SKDMrOJqVQqISl9NSHNRbo4/dv00rFSqVTrUK0KReRvPEhqjoiBAuIpRHt7e6xatarWYZhZanAQurq2snHjBk48sYe77noPra1tLF06iaamWkdnQyStjoj2yv1VHwBMK74aWAi8g+QW1U+AS4H/HOsgzawxdHfD5s19rF3bQUvLAP39l9DRsZ7u7oOYO7fW0dlostyquhHYApwMvDd9/+0igzKziW3NGpg9+1ZaWpIbFy0tA8yZs5Te3trGZdlkSRyviojLIuKX6esfgH0KjsvMJrBZs2DFipPo709uevT3N9PT08XMmbWNy7LJkjh+IOkUSbulr/cD2z0caGaWVWcntLa2MWPGvSxYcDkdHetpa2uls7PWkVkWo3aOS3qGZETVi+mu3YCt6fuIiMnFhTc23DluVn8GB2HKlNPp7z+MJUs+SWcn7hivMyN1jo/a4oiIvSJit4hoTl+7pfv22hWShpnVj/LhuM3N4sknb+DZZy9g3rxk28Nxdw2ZhuNKOgKYRtkorIi4pbiwxpZbHGZm+e1wi0PSYmAxyaiqeekr04A5ScdL2iBpk6QLhjkuSVenxx+UdGSWupL+Oj22XtKVWWIxM7OxMepzHMDREXFo3hNLagK+ChwH9AH3SVoWET8vK9YJTE9fHcA1QEe1upL+HJgPHBERf0indzczs3GSZVTVzyTlThzAUcCmiHgkIl4geR5kfkWZ+cD1kVgJ7CNp/1Hqng1cERF/AEindzczs3GSJXFcR5I8NqS3k9ZKyrLqXyvwaNl2X7ovS5lqdQ8CjpV0j6S7Jf3xcB8u6UxJqySt2rJlS4Zwzcwsiyy3qhaTTKO+lm1DcrPQMPsqe+JHKlOtbjOwL3A08MfATZLeGBW9/BGxCFgESed4jrjNzKyKLInjNxGxbAfO3QccULbdBjyWsczuVer2AbekieJeSS8C+5FMhWJmZgXLcqvq4XQVwFMlvWfolaHefcB0SQdK2h04BahMQMuA09PRVUcDT0XE46PUvRV4J4Ckg0iSzO8yxGNmZmMgS4vjFcAfgNll+wKo+hxHRAxIOgfoAZqAxRGxXtJZ6fFrSVYUPAHYBDwHnFGtbnrqxcBiSeuAF4APV96mMjOz4uzQehy7Gj8AaGaW386sx7EH8FHgMGCPof0R8ZExjdDMzHYJWfo4bgBeC8wB7ibpqH6myKDMzKx+ZUkcb4qIi4GtEXEd8BfA4cWGZWZm9SpL4uhP//6XpBnA3iQTHpqZWQPKMqpqkaR9gU+RDIndE7i40KjMzKxujZo4IuLr6dsfAW8sNhwzM6t3I96qkvTfJVU7/keS3lFMWGZmVq+qtTheDayRtBpYTTKlxx7Am4A/JXlae7s1NszMbGIbMXFExJckfYVkeo+3A0cAzwMPAadFxG/GJ0QzM6snVfs4ImIQuCN9mZmZZRqOa2Zm9hInDjMzy2XUxCHpwCz7zMysMWRpcdw8zL7/O9aBmJnZrmHEznFJh5DMiLt3xcJNkymbJdfMzBpLtVFVBwNzgX2AeWX7nwE+VmBMZmZWx6o9x/Fd4LuSjomIn41jTGZmVseyTHK4SdJFJDPivlTeCzmZmTWmLJ3j3yWZSv1O4HtlLzMbI6VSCUmjvkqlUq1DNRt9zXFJvRExc3zCKYbXHLddyeAgTJlyOv39M1iyZAGdndDUVOuorBGNtOZ4lhbHbZJOKCAmM6swOAhdXVuZOvVczj57kIULN9DVtZXBwVpHZrZNlj6OvwUukvQC8AIgICJicqGRmTWg7m7YvLmPtWs7aGkZoL//Ejo61tPdfRBz59Y6OrPEqC2OiNgrInaLiD0iYnK67aRhVoA1a2D27FtpaRkAoKVlgDlzltLbW9u4zMplmXJE6aJOF6fbB0g6qvjQzBrPrFmwYsVJ9PcnNwP6+5vp6eli5szaxmVWLksfx9eAY4APptvPAl8tLCKzBtbZCa2tbcyYcS8LFlxOR8d62tpa6eysdWRm22RJHB0R8VfA7wEi4klg90KjMmswQ8Nxm5vF8uWT2bjxYq66aitr1nyc5csn09zs4bhWP7IMx70HeBtwX0QcKWkKsCIiZo1HgGPBw3HNzPLbmeG4VwNLgamS/hH4CXD5GMdnZma7iKrDcSXtBvwSWAC8i2Qo7kkR8dA4xGZmZnVotDXHX5T0+Yg4Bnh4nGIyM7M6luVW1QpJJ0tS4dGYmVndy/Lk+MeBScCApN/jJ8fNzBpalj6O4yPi38YpHjMzq3NVb1VFxIvA58YpFjMz2wUU2sch6XhJGyRtknTBMMcl6er0+IOSjsxR9+8khaT98sZlZmY7rrA+DklNJFOTHAf0AfdJWhYRPy8r1glMT18dwDVAx2h1JR2QHvtN5m9qZmZjIs/suLvnnB33KGBTRDwSES8ANwLzK8rMB66PxEpgH0n7Z6j7BZJnS6o/9m5mZmNu1BaHpD8Zbn9E/GiUqq3Ao2XbfSStitHKtFarK+lEYHNEPFDt7pmkM4EzAV7/+tePEqqZmWWV5VbVJ8re70HSGlgNvHOUesP9X72yhTBSmWH3S3ol8PfA7FE+m4hYBCyCZK6q0cqbmVk2oyaOiJhXvp32L1yZ4dx9wAFl223AYxnL7D7C/j8CDgSGWhttwP2SjoqI32aIyczMdlKWUVWV+oAZGcrdB0yXdKCk3YFTgGUVZZYBp6ejq44GnoqIx0eqGxFrI2JqREyLiGlpLEc6aZiZjZ8sfRxfZtstpt2AmcADo9WLiAFJ5wA9QBOwOCLWSzorPX4tcDtwArAJeA44o1rdfF/NzMyKkGU9jg+XbQ4Av9rVniT3ehxmZvmNtB7HiC2OdMGmKRFxXcX+wyRNiYgtBcRpZmZ1rlofx5eBKcPsbwO+VEw4ZmZW76oljsMj4u7KnRHRAxxRXEhmZlbPqiWOlh08ZmZmE1i1xPELSSdU7pTUCTxSXEhmZlbPqg3HPQ+4TdL7SZ4UB2gHjgHmFh2YmZnVpxFbHBGxETgcuBuYlr7uBo5Ij5mZWQOq+gBgRPwB+JdxisXMzHYBOzLliJmZNTAnDjMzyyVz4pDUImmWpKlFBmRmZvVtxMQh6VpJh6Xv9yaZ2PB6YI2kU8cpPjMzqzPVWhzHls1IewawMSIOB95KsmyrmZk1oGqJ44Wy98cBtwJ47QsbD6VSCUnpqwlpLtLF6d+ml46VSqVah2rWcEacVl3SD4DPA5uBHwCHRMRvJTUD6yLikPELc+d4WvVd1+AgdHVtZePGDZx4Yg933fUeWlvbWLp0Ek1NtY7ObGLLPa068D+Aq4HXAueWtTTeBXxv7EM02153N2ze3MfatR20tAzQ338JHR3r6e4+iLmev8CsJqo+OR4Rx0fEzIj4Rtn+nog4f1yis4a3Zg3Mnn0rLS0DALS0DDBnzlJ6e2sbl1kjq7aQ09XVKkbE34x9OGYvN2sWLFx4Epde+qm0xdFMT08Xl15a68jMGle1W1VnAeuAm4DHAI1LRGZlOjth0aI2Zsy4l/nzv8+dd55MW1srnZ21jsyscVUbVbU/sAiYA5xGsgbHsoi4rnI5WbOxNjSqqrlZLF8+mY0bL+aqq7ayZs3HWb58Ms3NHlVlVisjjqp6WSGpFTgV+DjwyYi4oejAxpJHVZmZ5bcjo6qGKh5JkjSOA7rZtjaHmZk1oGqd458mWbDpIeBG4MKIGBivwMzMrD5Va3FcTLJE7FvS1+WSIOkkj4g4ovjwzMys3lRLHAeOWxRmZrbLGDFxRMSvh9svqQk4BRj2uJmZTWzVplWfLOlCSV+RNFuJvya5ffX+8QvRzMzqSbVbVTcATwI/A/4S+ASwOzA/InqLD83MzOpRtcTxxnT9DSR9Hfgd8PqIeGZcIjMzs7pU7cnx/qE3ETEI/NJJw8zMqrU43iLp6fS9gFek20PDcScXHp2ZmdWdaqOqvEyOmZltp9qtKjMzs+04cZiZWS6FJg5Jx0vaIGmTpAuGOS5JV6fHH0wnVKxaV9JVkh5Oyy+VtE+R38HMzF6usMSRPmH+VaATOBQ4VdKhFcU6genp60zgmgx17wBmpHNlbQQuLOo7mJnZ9opscRwFbIqIRyLiBZIZdudXlJkPXB+JlcA+kvavVjciVpTN0rsSaCvwOzSkwUG47Ta47LLk7+BgrSMys3pSZOJoBR4t2+5L92Upk6UuwEdI1gjZjqQzJa2StGrLli05Q288QyvuSU1MnryM889fzTPPfIbzz1/N5MnLkJq84p6ZAcUmjuHWKK9cbnCkMqPWlfT3wADwzeE+PCIWRUR7RLRPmTIlQ7iNrVQqEREsXz7IIYcczLp1R3PllRexbt3RHHzwISxfPkhEOHGYWaGJow84oGy7DXgsY5mqdSV9mGSRqQ9FlrVvLbM1a2D27FtpaUnuBra0DDBnzlJ6e2sbl5nVjyITx33AdEkHStqdZCr2ZRVllgGnp6OrjgaeiojHq9WVdDzwSeDEiHiuwPgb0qxZsGLFSfT3J8+G9vc309PTxcyZtY3LzOpHYYkj7cA+B+ghWX72pohYL+ksSWelxW4nmaZ9E/BPwP+sVjet8xVgL+AOSb2Sri3qOzSSoT6OefOaePjhDRx22D0sWHA5hx12Dxs2PMy8ee7jMLOEGuFOT3t7e6xatarWYewyBgehuxt6e2HmTOjshCZPQGPWcCStjoj2yv3VJjm0BtXUBHPnJi8zs0qecsTMzHJx4jAzs1ycOMzMLBcnDjMzy8WJw8zMcnHiMDOzXJw4zMwsFycOMzPLxYnDzMxyceIwM7NcnDjMzCwXJw4zM8vFicPMzHJx4jAzs1ycOMzMLBcnDjMzy8WJw8zMcnHiMDOzXJw4zMwsFycOMzPLxYnDzMxyceIwM7NcnDjMzCyX5loHMNEMDkJ3N6xZA7NmQWcnNDXVOiozs7HjFscYGhyErq6tLFy4geee+ywLF26gq2srg4O1jszMbOw4cYyBUqmEJJqb57Jhw8OsXDmDz3zmAlaunMHDD2+guXkukiiVSrUO1cxsp/lW1RgolUqUSiUuuwyee+6ztLQMANDSMsDJJ9/BpEm38alP1ThIM7Mx4hbHGBhqcVxyyVxuueXd9Pcn+bi/v5mbbz6Oiy92i8PMJg5FRK1jKFx7e3usWrWq8M8Z6uPo69vMnDlL6enpoq2tlaVLJ7mD3Mx2OZJWR0T7dvudOMbW0Kiq3l6YOdOjqsxs1zVS4nAfxxhraoK5c5OXmdlE5D4OMzPLxS2OEfhBPjOz4RXa4pB0vKQNkjZJumCY45J0dXr8QUlHjlZX0qsk3SHpF+nffcc6bj/IZ2Y2ssISh6Qm4KtAJ3AocKqkQyuKdQLT09eZwDUZ6l4A3BUR04G70u0x4Qf5zMxGV2SL4yhgU0Q8EhEvADcC8yvKzAeuj8RKYB9J+49Sdz5wXfr+OuCksQ99FvPnr3jZg3wnndQDzBz7jzIz28UUmThagUfLtvvSfVnKVKv7moh4HCD9O3W4D5d0pqRVklZt2bIlU8ClUomIYPnyy7jrrve87EG+O+88meXL/4GIcIvDzBpakYlDw+yrfGhkpDJZ6lYVEYsioj0i2qdMmZKnKp2d0NraRkfHei688Ao6OtbT1tZKZ2eu05iZTUhFjqrqAw4o224DHstYZvcqdf9D0v4R8Xh6W+uJMY2aZPTU0qWT6O4+iN7eT3LppR5VZWY2pMjEcR8wXdKBwGbgFOCDFWWWAedIuhHoAJ5KE8KWKnWXAR8Grkj/freI4P0gn5nZ8ApLHBExIOkcoAdoAhZHxHpJZ6XHrwVuB04ANgHPAWdUq5ue+grgJkkfBX4DvK+o72BmZtvzXFVmZjaskeaq8pQjZmaWixOHmZnl4sRhZma5NEQfRzpK69c5q+0H/K6AcHZWvcYF9RtbvcYFjm1H1GtcMPFie0NEbPcgXEMkjh0hadVwnUK1Vq9xQf3GVq9xgWPbEfUaFzRObL5VZWZmuThxmJlZLk4cI1tU6wBGUK9xQf3GVq9xgWPbEfUaFzRIbO7jMDOzXNziMDOzXJw4zMwsl4ZLHEWsgz6OsX0ojelBST+V9JayY7+StFZSr6QxnZgrQ1x/Jump9LN7JV2Ste44xPaJsrjWSRqU9Kr0WJG/2WJJT0haN8LxWl5no8VWq+tstLhqeZ2NFlutrrMDJP1A0kOS1kv622HKjP21FhEN8yKZafffgTeSrPnxAHBoRZkTgG6SxaSOBu7JWnccYnsbsG/6vnMotnT7V8B+NfrN/gy4bUfqFh1bRfl5wL8W/Zul5/4T4Ehg3QjHa3KdZYxt3K+zjHHV5DrLElsNr7P9gSPT93sBG8fj/2mN1uIoah30cYktIn4aEU+mmytJFrgq2s5875r/ZhVOBZaM4eePKCJ+BPy/KkVqdZ2NGluNrrMsv9lIav6bVRjP6+zxiLg/ff8M8BDbL9E95tdaoyWOotZBH6/Yyn2U5F8RQwJYIWm1pDNrENcxkh6Q1C3psJx1i44NSa8EjgduLttd1G+WRa2us7zG6zrLqhbXWWa1vM4kTQNmAfdUHBrza63IFQDrUU3XQR9F5vNL+nOS/6DfUbb77RHxmKSpwB2SHk7/lTQecd1PMqfNs5JOAG4FpmesW3RsQ+YB/xYR5f9qLOo3y6JW11lm43ydZVGr6yyPmlxnkvYkSVbnRsTTlYeHqbJT11qjtTh2Zh30LHWLjg1JRwBfB+ZHxH8O7Y+Ix9K/TwBLSZqh4xJXRDwdEc+m728HWiTtl6Vu0bGVOYWK2wcF/mZZ1Oo6y6QG19moanid5THu15mkFpKk8c2IuGWYImN/rRXRYVOvL5IW1iPAgWzrDDqsosxf8PKOpHuz1h2H2F5Psszu2yr2TwL2Knv/U+D4cYzrtWx7mPQokiV9VQ+/WVpub5L705PG4zcr+4xpjNzRW5PrLGNs436dZYyrJtdZlthqdZ2l3/964ItVyoz5tdZQt6qiuHXQxyu2S4BXA1+TBDAQyWyXrwGWpvuagW9FxPfHMa73AmdLGgCeB06J5Mqsh98MoAtYERFby6oX9psBSFpCMgpoP0l9wEKgpSyumlxnGWMb9+ssY1w1uc4yxgY1uM6AtwOnAWsl9ab7LiJJ/oVda55yxMzMcmm0Pg4zM9tJThxmZpaLE4eZmeXixGFmZrk4cZiZWS5OHLbLS2ciHZqVdLmkfdL90ySFpMvKyu4nqV/SV9LtgyX9MK3/kKTtVkmT9EtJB1fs+6KkBVVierYshg+O0VetStK5kk7fwbo/lNSevr8oQ/kbJU3fkc+yXZ8Th00Ez0fEzIiYQfIA1l+VHXsEmFu2/T6gfKz61cAX0vpvBr48zPlvJHkiGABJu5E8U/DtDLFNAwpPHJKagY8A3xqD042aOIBrgBETp01sThw20fyMl0/U9jzw0NC/poEPADeVHd+fZOoFACJi7TDnXEJZ4iCZYvtXEfFrSR9PWzrrJJ07TN0rgGPTFs15aQvkx5LuT19vgyQZSfqakjUVbpN0u6T3psfeKunudJK8nnRm00rvBO5PH+p6s6R7hw6kn/lg+v5dktYoWR9isaT/Vn4SSVcAr0jj/aakSZK+l04suE7SB9KiPwbenSYsazBOHDZhSGoC3gUsqzh0I3CKpDZgkJfPx/MF4F+VzLZ63tBtrnIR8SDworYtaHQKsETSW0mewu0gmcrhY5JmVVS/APhx2qL5AvAEcFxEHEmSxK5Oy72HpHVyOPCXwDHpd2ohaQW9NyLeCiwG/nGYr/92YHUa70PA7pLemB77AHCTpD2AbwAfiIjDSZ5kPrviu17Athbch0hmen0sIt6Stui+n5Z7keRJ5LdgDceJwyaCV6TTLfwn8Crgjorj3weOI1kn4WW3lyLiX4A3A98hmVJiZeW/wlNLSJJPM8maBd8hmTV2aURsjWTyvVuAY0eJtQX4J0lr03Mcmu5/B/CdiHgxIn4L/CDdfzAwg2RW1V7gUwy/Psb+wJay7ZuA96fvP5B+74OBX0bExnT/dSStp2rWkrQsPivp2Ih4quzYE8DrRqlvE5ATh00Ez0fETOANJJO1lfdxEMkiNauB83n5OglDxx+LiMURMR8YIPkfdaUlJP8jfjfwYCQznQ43LfVozgP+g+Rf6u1pvFQ5l4D1aQtgZkQcHhGzhyn3PLBH2fa3gfdLOgiIiPjFjsSbJpm3kiSQz6hsudb0857Pe07b9Tlx2ISR/mv4b4C/S2/xlPs88MkomyIcXlpzuSV9/1qSyf02D3Pufydp0VzBtmmzfwScJOmVkiaRTHL344qqz5As6Tlkb+Dx9FbPaSSTywH8BDg57et4DUnrB2ADMEXSS7eutG0Bo3IPAW+qiHcQuJhtrayHgWmShsqdBtw9zLn6y36T1wHPRcT/AT5HsnzqkIN4+UADaxDu2LIJJSLWSHqApB/ix2X71zP8/+RmA1+S9Pt0+xPpraLhLAE+Q7KmAhFxv6RvAEMd0V+PiDUVdR4kmb31AZL+ha8BN0t6H8ntqKGZVG8m6Z9ZR7Ju9D3AUxHxQtpJfrWkvUn+m/3iMN+lG7ihYt+3gatIps0mIn4v6QzgO+ktt/uAa9neIuBBSfeTTNl9laQXgX7SPpE0uT0fEY+P8FvZBObZcc3qhKQ9I1nd7tUkyejtVZLYcPWXAgvS21KFknQe8HRE/HPRn2X1xy0Os/pxWzqqa3fgsjxJI3UBSSd54YkD+C+2b+FYg3CLw8zMcnHnuJmZ5eLEYWZmuThxmJlZLk4cZmaWixOHmZnl8v8B7wUjUcLw1igAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "plt.errorbar(Vrms, Irms, errIrms, fmt = 'ko', markersize = 5,\\\n", " linewidth = 1.8,\\\n", " markeredgecolor = 'b',\\\n", " markerfacecolor = 'yellow',\\\n", " capsize = 5)\n", "plt.xlabel('RMS Voltage (volts)')\n", "plt.ylabel('RMS Current (amps)');" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "To do the actual fit, we will use the 'curve_fit()' function from, the *SciPy* module. This way of fitting is very nice because we will be able to use it for all types of fit models (linear, polynomial, linear-in-parameter fits, and nonlinear fits). As we will see, it is capable of doing both unweighted and weighted fits and it will return uncertainties in the fit parameters." ] }, { "cell_type": "code", "execution_count": 77, "metadata": {}, "outputs": [], "source": [ "from scipy.optimize import curve_fit" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The first step is to define a function for the model that we will fit our data to. In this case, the model is just the equation of a staight a line." ] }, { "cell_type": "code", "execution_count": 78, "metadata": {}, "outputs": [], "source": [ "def linearFunc(x, intercept, slope):\n", " y = slope*x + intercept\n", " return y" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here is the actual command to execute the fit. At a minimum, *curve_fit()* requires as inputs the function that defines the model, the $x$-data, and the $y$-data. The statement below tells *curve_fit()* to return a list of the the best-fit parameters (a_fit) and the covariance matrix (cov) which will be used to determine the error in the fit parameters." ] }, { "cell_type": "code", "execution_count": 79, "metadata": {}, "outputs": [], "source": [ "a_fit, cov = curve_fit(linearFunc, Vrms, Irms)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here are the contents of a_fit and cov." ] }, { "cell_type": "code", "execution_count": 80, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The best-fit parameters are: (1) slope = 0.004562609977167653 and (2) intercept -3.029377062510313e-05\n", "Here is the covariance matrix:\n", " [[ 2.10867163e-10 -1.45419351e-10]\n", " [-1.45419351e-10 1.57145259e-10]]\n" ] } ], "source": [ "print('The best-fit parameters are: (1) slope =', a_fit[1], 'and (2) intercept',\\\n", " a_fit[0])\n", "print('Here is the covariance matrix:\\n', cov)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As we will see in PHYS 232, the uncertainties of the best-fit parameters are determined from the square roots of the diagonal elements of the covariance matrix. We can select the diagonal elements using:" ] }, { "cell_type": "code", "execution_count": 81, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[2.10867163e-10 1.57145259e-10]\n", "The error in the slope is: 1.2535759226262432e-05 \n", "The error in the intercept is: 1.4521265884900987e-05\n" ] } ], "source": [ "print(np.diag(cov))\n", "print('The error in the slope is:', np.sqrt(np.diag(cov))[1], '\\n'\n", " 'The error in the intercept is:', np.sqrt(np.diag(cov))[0])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "*NumPy* has a nice package for polynomials, called *polynomial*. There are six different polynomial types in this package. For our case, we are dealing with a simple power series. You can use the*Polynomial* constructor for this. y = np.polynomial.Polynomial([a, b, c]) results in $y = a + b\\,x + c\\,x^2$. \n", "\n", "Here's the best-fit fucntion obtained using *a_fit* from *curve_fit()* and the built in polynomial package of *NumPy*." ] }, { "cell_type": "code", "execution_count": 82, "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$x \\mapsto \\text{-3.029377062510313e-05} + \\text{0.004562609977167653}\\,x$" ], "text/plain": [ "Polynomial([-3.02937706e-05, 4.56260998e-03], domain=[-1, 1], window=[-1, 1])" ] }, "execution_count": 82, "metadata": {}, "output_type": "execute_result" } ], "source": [ "fitFcn = np.polynomial.Polynomial(a_fit)\n", "fitFcn" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can then evaluate *fitFcn(x)* at any value of $x$ that we want." ] }, { "cell_type": "code", "execution_count": 83, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(-3.029377062510313e-05, 0.004532316206542549)" ] }, "execution_count": 83, "metadata": {}, "output_type": "execute_result" } ], "source": [ "fitFcn(0), fitFcn(1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This gives us an easy way to plot the fit on top of the data." ] }, { "cell_type": "code", "execution_count": 84, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEWCAYAAABxMXBSAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA4QUlEQVR4nO3deZxN9f/A8dc7pk2lksqXikplZjCWDGnRQijkV99vtOk7Q1LIN0tKMRUhEb52RiVFdsIY8rUmMtYxtrTZo5I1meX9++Ocqds0c+eO5s6Z5f18PO7Dved8Pue+7+007/v5fM75fERVMcYYYwJ1jtcBGGOMKVgscRhjjMkRSxzGGGNyxBKHMcaYHLHEYYwxJkcscRhjjMkRSxymyBORx0VkYYBlnxaRlUGM5ayOLyKviMi4YMSUW0QkSUTqeR2H+fsscRhPiYiKyI0ZtsWIyMS8ikFVP1LVBrlxLBFZKiKtc+NYmRy7vPt9Fc+4T1XfUtWgvG9Oud/BaRE54fOoo6phqrrULZOn/41N7rLEYYw5ayJSLItd7VX1Ip/HF3kamAkqSxwmXxOReiKyV0Q6i8ghETkgIv9291UQkV9E5Bz39TgROeRTd6KIdHKflxSRWLf+PhHpnf5HL2P3kIg0EJEdInJUREaIyLKMrQgReUdEjojItyLSyN3WB7gDGOb+yh7mbr9FRBaJyM/ucf/lc5xSIjJHRI6JyJfADWf5Pf3+C96nZdJKRHaLyI8i0sOn7Dki0l1EvhaRn0Rkiohc7rN/qogcdD//chEJ89n3voiMFJH5InISuDsHMX4nIveJSEPgFeBR93vadDaf2XjHEocpCK4GSgJlgWhguIhcpqrfAseAam65O4ATIlLJfX0nsMx9/gGQAtzolm8A/KVrR0SuAKYBLwOlgB3AbRmKRbrbrwDeBmJFRFS1B7CCP35ttxeREsAi4GPgSqAlMMLnj/Fw4DRQBohyH7nlduBm4F6gp8/30hF4CLgL+AdwxI0jXRxQ0Y13PfBRhuM+BvQBLgZyPB6jqguAt4BP3O+pak6PYbxlicMUBMnAG6qarKrzgRM4fxDBSQx3icjV7utp7usKwCXAJhG5CmgEdFLVk6p6CHgXaJHJezUGklR1hqqmAEOBgxnKfK+qY1U1FSchlQGuyiL2B4HvVPU9VU1R1fXAdOARt8XzMNDTjWuLe7zc8rqq/qqqm4BNQPof6LZAD1Xdq6q/ATFuPMUBVHW8qh732VdVREr6HHe2qn6uqmmqejqL9x7qtgZ/EZH1ufiZTD7wl0E2Y/JYKhCSYVsITrJI95P7RzzdKeAi9/kyoCmwF1gOLAWexPkVv0JV00TkOveYB0Qk/RjnAHsyiecfvttVVUVkb4YyB332n3KPeRGZuw6IFJFffLYVBz4ESrvPfeP4PovjnA3fhOf7nV0HzBSRNJ/9qcBVInIQpzXxTze+9DJXAEfd55l9bxl1VNV8fZWXOXuWOIzXdgPlgW0+2yoAOwOsvwwYgJM4luF0nYzCSRzp3VR7gN+AKzIkoMwcAMqlvxAnK5TLuvhfZJxueg+wTFXrZyzotjhSgGuA7e7ma3PwXmdrDxClqp9nEtOTQDPgPuA7nC7CI4D4FMuNKbVtWu4CzLqqjNc+AV4VkXLuoO19QBOcLqdsqepXwK/AE8ByVT0G/IDTBbTMLXMAWAgMFJFL3Pe5QUTuyuSQ84DKIvKQ23XzPM4YS6B+AK73eT0XuElEnhSREPdxq4hUcru6ZgAxInKhiIQCrQJ4j/NE5HyfR07/Px4F9HFbYohIaRFp5u67GCfJ/gRciDMWEQw/AOXPInaTD9h/NOO1N4BVOC2FIziDzY+7/f2BWobTnbXb57UAG3zKPAWcC2x132caztjEn6jqjzjdNG/j/PEMBRJw/pgGYgjOeMERERmqqsdxBuJbAPtxuo/6A+e55dvjdCEdBN4H3gvgPU7gJMv0xz0BxuYb4xxgoYgcB1bjDPgDTMDpLtuH812tzuGxAzXV/fcnGwMpeMQWcjIma+4v4r04yWyJ1/EYkx9Yi8OYDETkfhG5VETOw7nfQAjeL29jCpygJg4Raeje8LRLRLpnsl9EZKi7f7OIVPfZN16cG762ZKhzuXsz1Vfuv5cF8zOYIqkO8DXwI854y0Oq+qu3IRmTfwStq8q9YmQnUB+nqb8WaKmqW33KNAY64Fw7HwkMUdVId9+dOH25E1Q13KfO28DPqtrPTUaXqepLQfkQxhhj/iKYLY5awC5V/UZVzwCTcS7z89UMJzGoqq4GLhWRMgCquhz4OZPjNuOPm6Q+wLkD1hhjTB4J5n0cZfnzjUJ7+ePKDX9lyuJcS5+Vq9zLK1HVAyJyZWaFROQZ4BmAEiVK1LjllltyFr0xxhRBJ06c4LvvvuO3334D+FFVS2csE8zEIZlsy9gvFkiZs6KqY4AxADVr1tSEhITcOKwxxhRKx48fp3v37owYMYIKFSowbtw47r333kxnMghmV9VenDti05XDuY49p2Uy+iG9O8v991A25Y0xxvixYMECwsLCGDlyJJ06dSIxMZF77sn69qBgJo61QEVxpr4+F+cGqDkZyswBnnKvrqoNHE3vhvJjDn/cXdsKmJ2bQRtjTFHx888/8/TTT9OoUSMuuugiPv/8c959911KlCjht17QEoc7J1B7IB5nHqIpqpokIs+KyLNusfnAN8AuYCzwXHp9EZkEfAHcLM56DNHurn5AfRH5CueKrX7B+gzGGFNYTZ8+ndDQUD766CNeffVVNmzYQJ06dQKqWyTuHLcxDmOMcRw8eJD27dszffp0qlevTmxsLBEREZmWFZF1qloz43a7c9wYY4oAVWXChAmEhoYyd+5c+vbty5o1a7JMGv7YtOrGGFPI7d69m7Zt27JgwQLq1q3LuHHj+Du3KFiLwxhjCqm0tDRGjhxJWFgYK1asYOjQoSxfvvxvJQ2wFocxxhRKX331Fa1bt2b58uXUr1+fMWPGUL58+Vw5trU4jDGmEElJSWHAgAFUqVKFTZs2MX78eOLj43MtaYC1OIwxptBITEwkKiqKhIQEHnroIYYPH84//vGPXH8fa3EYY0wBd+bMGWJiYqhRowbff/89n3zyCTNmzAhK0gBrcRhjTIG2du1aoqKi2LJlC48//jiDBw/miiuuCOp7WovDGGMKoFOnTtGlSxdq167NkSNHmDt3LhMnTgx60gBrcRhjTIGzbNkyWrduza5du2jbti39+/enZMmSefb+1uIwxhiPxcTEICLZPl5++WXatWtHvXr1SEtL43//+x+jRo3K06QBNleVMcbkO+mD2vv3/7HKxPz582nbti379++nU6dOvPnmm1x44YVBjcPmqjLGmAIgNRVOn76P48c7MXcuHDr0E08++SQPPPAAl1xyCatWrWLgwIFBTxr+2BiHMcbkE6mp0Lz5Sa68shNNm8bTo0cSu3fv4fjxKbz22mv06NGD8847z+swLXEYY0x+ERcH+/btJTExkpCQFPr06Unlymt5441tdOhwvdfh/c66qowxJp9Yv16pX38WISEpAISEpPDQQ/EcPZp/kgZY4jDGmHzh+++/Z9asXsyceR/JyU5nUHJycRYubM5ZLJkRVJY4jDHGQ2lpaQwfPpywsDB27BhMiRIXER7+Jd26vUVkZBLlypWlUSOvo/wzSxzGGOORnTt3ctddd9G+fXvq1q3L1q2JrF17M4cPv8vIkefwxhs3MXNmCYoV8zrSP7PEYYwxeSwlJYW3336bKlWqsGXLFpo1a8bChQspX748xYsLR458yIkT3WnSRChe/I8bAGNiYrwOHbCrqowxJk9t3ryZqKgo1q1bR/PmzRk+fDhlypTxOqwcsRaHMcbkgd9++42ePXtSo0YN9uzZw5QpU5g+fXqBSxpgLQ5jjAm6NWvWEBUVxdatW3niiScYPHgwpUqV8jqss2YtDmOMCZJTp07RuXNnbrvtNo4dO8a8efP48MMPC3TSAGtxGGNMUCxdupTWrVvz9ddf07ZtW95++20uueQSr8PKFdbiMMaYXHTs2DGeffZZ7r77bgCWLFnCqFGjCk3SAEscxhiTa+bPn09YWBhjx46lc+fObN68mXr16nkdVq6zxGGMMX/TTz/9xBNPPMEDDzxAyZIl+eKLL3jnnXc8nfo8mCxxGGPMWVJVpkyZQqVKlfjkk0/o2bMn69ato1atWl6HFlQ2OG6MMWfhwIEDPPfcc8yaNYsaNWrw2WefUaVKFa/DyhPW4jDGmBxQVd577z1CQ0NZsGABb7/9NqtXry4ySQOsxWGMMQH77rvvaNu2LQsXLuSOO+5g3Lhx3HTTTV6HleesxWGMMdlIS0tj2LBhhIeHs2rVKoYPH87SpUuLZNIAa3EYY4xfO3bsoHXr1qxcuZL777+f0aNHc91113kdlqesxWGMMZlISUmhf//+VK1alaSkJN5//33i4uKKfNKAICcOEWkoIjtEZJeIdM9kv4jIUHf/ZhGpnl1dEYkQkdUislFEEkSkcF/3ZozJdTExMb+vcZHVIyQkhO7du/PAAw+wdetWWrVqhYh4HXq+ELTEISLFgOFAIyAUaCkioRmKNQIquo9ngJEB1H0beF1VI4Ce7mtjjAlYTEwMqoqqkpKiXHbZk1x0UX9mzDjDK6+8RvHixbnqqquYNm0a06dP5+qrr/Y65HwlmGMctYBdqvoNgIhMBpoBW33KNAMmqKoCq0XkUhEpA5T3U1eB9ElfSgL7g/gZjDGFWGoqNG9+kiuv7ETTpvG88soWdu+uSYsWjzNkyCAuv/xyr0PMl4KZOMoCe3xe7wUiAyhTNpu6nYB4EXkHp8V0W2ZvLiLP4LRiuPbaa8/qAxhjCre4ONi3by+JiZGEhKTQp09PatbcxKOPvo/ljKwFc4wjs85ADbCMv7rtgP+o6jXAf4DYzN5cVceoak1VrVm6dOkAQzbGFCUzZnzDvffOICQkBYCQkBQaN/6UjRu9jSu/C2bi2Atc4/O6HH/tVsqqjL+6rYAZ7vOpOF1ixhgTsKNHj9K2bVvee68jc+bcT3Ky0/mSnFyc+PjmRER4G19+F8zEsRaoKCIVRORcoAUwJ0OZOcBT7tVVtYGjqnogm7r7gbvc5/cAXwXxMxhjCpm5c+cSFhbGuHHj6Nw5jBtvrEh4+Jd06/YWkZFJlCtXlkaNvI4yfwvaGIeqpohIeyAeKAaMV9UkEXnW3T8KmA80BnYBp4B/+6vrHroNMEREigOncccxjDHGnx9//JEXXniBjz/+mPDwcGbOnMmtt95KaiqULv0uI0eGMWnSTTRqBMWKeR1t/ibOBU2FW82aNTUhIcHrMIwxHkif+rxDhw788ssv9OjRg5SUFHr37p1t3V69ehETExP8IPMpEVmnqjUzbrcpR4wxhdb+/ft57rnnmD17NrfeeiuxsbFUrlwZgDfffNPj6Aoum3LEGFPoqCqxsbGEhoYSHx/PgAEDWLVq1e9Jw/w91uIwxhQq3377Lc888wyfffYZd955J+PGjaNixYpeh1WoWIvDGFMopKamMnToUMLDw1m9ejUjR45kyZIlljSCwFocxpgCb/v27URHR7Nq1SoaNmzI6NGjbcaIILIWhzGmwEpOTqZv375ERESwfft2JkyYwPz58y1pBJm1OIwxBdKGDRuIjo5mw4YNPPLIIwwbNoyrrrrK67CKBGtxGGMKlNOnT9OjRw9uvfVW9u/fz/Tp05k6daoljTxkLQ5jTIGxatUqoqOj2b59O08//TQDBw60qc89YC0OY0y+d/LkSTp16sTtt9/OqVOnWLBgAe+9954lDY9Yi8MYk6999tlntGnThu+++47nnnuOfv36cfHFF3sdVpHmN3GIyPnAg8AdwD+AX4EtwDyfSQeNMSbX/fLLL3Tp0oXY2FgqVqzI8uXLueOOO7wOy+Cnq0pEYoDPgTrAGmA0MAVIAfqJyCIRqZIXQRpjipY5c+YQFhbGe++9R7du3di0aZMljXzEX4tjrarGZLFvkIhcCdjF0saYXHP48GE6duzI5MmTqVy5MrNnz6Zmzb9Mzmo8lmWLQ1XnZdwmIueIyCXu/kOqanOVG2P+NlVl0qRJhIaGMn36dF5//XUSEhIsaeRT2V5VJSIfi8glIlIC2ArsEJGuwQ/NGFMU7Nu3j6ZNm/LYY49x/fXXs379enr27Mm5557rdWgmC4FcjhuqqseAh3BW7LsWeDKYQRljCj9VZezYsYSGhrJ48WLeeecdVq1aRXh4uNehmWwEcjluiIiE4CSOYaqaLCKFf9lAY0zQfPPNN7Rp04b//e9/1KtXj7Fjx3LjjTd6HZYJUCAtjtHAd0AJYLmIXAccC2ZQxpjCKTU1lcGDB1O5cmXWrl3LqFGjWLx4sSWNAibbFoeqDgWG+mz6XkTuDl5IxpjCaNu2bURHR/PFF1/QuHFjRo0axTXXXON1WOYsBDI4XkpEhorIehFZJyJDgJJ5EJsxphBITk6mT58+REREsGPHDj788EPmzp1rSaMAC6SrajJwGHgYeMR9/kkwgzLGFA4bNmygVq1avPrqqzRr1oytW7fyxBNPICJeh2b+hkASx+Wq+qaqfus+egOXBjkuY0wBdvr0aV555RVuvfVWDh48yIwZM5gyZYpNfV5IBJI4lohIC/fmv3NE5F/AX24ONMYYgM8//5yIiAj69u3LU089xdatW2nevLnXYZlcFEjiaAt8DJxxH5OBF0XkuIjY1VXGGABOnDhBx44dueOOOzh9+jTx8fGMHz+eyy67zOvQTC4L5Koqm7/YGANATEwMr7/+uvvqHKARUA3YAMQBaQC0bNmSBg0aeBKjCT5Rzf5ePncW3PL4JBpVnRG8sHJXzZo1NSHBptUyJrekpkLz5ifZuXMHTZvGM3t2A37++QjTpp3PXXfd7nV4JpeIyDpV/cuEYdm2OERkPFAFSCL95wQoUGAShzEmd8XFwffff0diYiQhISn06dOTyMgkjh+/yevQTB4IZMqR2qoaGvRIjDEFwqFDh3j55aU0avQ1ISEpAISEpHD//TPZuPElHnzQ4wBN0AUyOP6FiFjiMKaIU1U++ugjQkND2b59EnFxD5Kc7Pz2TE4uTnx8cyIivI3R5I1AEscHOMljh4hsFpFEEdkc7MCMMfnH3r17adKkCU888QQ33ngj69b1pkKF6wkP/5Ju3d4iMjKJcuXK0qiR15GavBBIV9V4nGnUE/ljjMMYUwSkT33etWtXkpOTGTRoEB07dqRYsWLMnAmlS7/LyJFhTJp0E40aQbFiXkds8kIgLY7dqjrHvWv8+/RH0CMzxnjq66+/5t5776Vt27bUqFGDxMREjh49SvHixRERihcXjhz5kBMnutOkifNaxHnExMR4Hb4JomwvxxWREThTjHwK/Ja+3S7HNaZwSk1NZejQofTo0YOQkBDeeecdWrdubfNLFUFZXY4bSIvjApyE0QBo4j4Cum5CRBq6YyO7RKR7JvvFnXl3lzt+Uj2QuiLSwd2XJCJvBxKLMSZ7W7du5fbbb+fFF1/k3nvvJSkpiTZt2ljSMH8SyJ3j/z6bA4tIMWA4UB/YC6wVkTmqutWnWCOgovuIBEYCkf7qumuBNAOqqOpvInLl2cRnjPlDcnIy/fr1o3fv3lx88cV89NFHtGzZ0hKGyVQgNwCeD0QDYcD56dtVNSqbqrWAXar6jXucyTh/8H0TRzNggjr9ZatF5FIRKYNzl3pWddsB/VT1NzeOQwF8TmNMFtatW0dUVBSbN2+mRYsWDBkyhCuvtN9jJmuBdFV9CFwN3A8sA8oBxwOoVxbY4/N6r7stkDL+6t4E3CEia0RkmYjcmtmbi8gzIpIgIgmHDx8OIFxjipZff/2V7t27ExkZyeHDh5k9ezaTJk2ypGGyFUjiuFFVXwNOquoHwANA5QDqZdbGzTgSn1UZf3WLA5cBtYGuwBTJpD2tqmNUtaaq1ixdunQA4RpTdKxcuZKIiAj69+/P008/zdatW2natKnXYZkCIpDEkez++4uIhOMsG1s+gHp7Ad+1IcsB+wMs46/uXmCGOr7EubfkigDiMabIO3HiBB06dODOO+/kzJkzLFq0iHHjxnHppZd6HZopQAJJHGNE5DLgVWAOzjhD/wDqrQUqikgFETkXaOHW9zUHeMq9uqo2cFRVD2RTdxZwD4CI3AScC/wYQDzGFGkLFy4kPDyc4cOH06FDBxITE7nvvvu8DssUQIFcVTXOfbocuD7QA6tqioi0B+KBYsB4VU0SkWfd/aOA+UBjYBdwCvi3v7ruoccD40VkC87CUq00kLnhjSmijhw5wosvvsj777/PzTffzMqVK7ntttu8DssUYFneACgiTwAfq2qm04yIyA1AGVVdGcT4coXdAGiKqpkzZ/Lcc89x+PBhXnrpJV577TXOP//87Csaw9mtx1EK2CAi64B1wGGcy3FvBO7C6R76y019xhjv/fDDD3To0IGpU6cSERHB/PnzqVatmtdhmUIiy8ShqkNEZBjOeEJdnMWcfgW2AU+q6u68CdEYE6j0qc9feOEFTpw4QZ8+fejatSshISFeh2YKEb9jHKqaCixyH8aYfGzPnj08++yzzJ8/nzp16hAbG0ulSpW8DssUQoFcVWWMycfS0tIYPXo0YWFhLF26lMGDB7NixQpLGiZoAlmPwxiTT+3atYs2bdqwdOlS7rnnHsaOHcv11wd88aMxZyXbFoeIVAhkmzEm76SmpjJw4ECqVKnC+vXrGTt2LJ999pklDZMnAumqmp7Jtmm5HYgxJjBbtmzhtttuo0uXLtx3331s3brV1ssweSrLrioRuQVnRtySIvJ/PrsuwWeWXGNM3jhz5szvU5+XLFmSSZMm8eijj1rCMHnO3xjHzTgLNl2Ks3hTuuNAmyDGZIzJICEhgaioKBITE2nZsiVDhgzBJu80XvF3H8dsYLaI1FHVL/IwJmOM69dff6VXr14MHDiQq6++mjlz5tCkSZPsKxoTRIFcVbVLRF7BmRH39/IBLORkjPkbVqxYQXR0NF999RWtW7dmwIABNoutyRcCSRyzgRXAZ0BqcMMxxhw/fpzu3bszYsQIKlSowGeffca9997rdVjG/C6QxHGhqr4U9EiMMcTHx/PMM8+wZ88eXnjhBfr06UOJEiW8DsuYPwnkcty5ItI46JEYU4T9/PPPPP300zRs2JALL7yQlStXMnjwYEsaJl8KJHG8gJM8TovIMRE5LiLHgh2YMUXF9OnTCQ0NZeLEifTo0YMNGzbYehkmXwtkIaeL8yIQY4qagwcP0r59e6ZPn061atVYsGABERERXodlTLYCmXJEROQJEXnNfX2NiNQKfmjGFE6qyoQJEwgNDWXu3Lm89dZbrFmzxpKGKTAC6aoaAdQBHnNfnwCGBy0iYwqx3bt307hxY1q1akWlSpXYuHEjL7/8sq2XYQqUQBJHpKo+D5wGUNUjwLlBjcqYQiYtLY2RI0cSHh7O8uXLGTJkCMuXL+eWW27xOjRjciyQy3GTRaQYoAAiUhrIdB1yY4qymJgYXn/9dffVOUAjoBqwAYgj/X+bjh070rFjR09iNCY3iKr6LyDyOPAoUB34AHgEeFVVpwY/vNxRs2ZNTUhI8DoMU0SkpkLz5ifZuXMHTZvGM3t2Aw4c+IGBAw/SuvW/bVJCU2CIyDpVrZlxu98Wh4icA3wLdAPuBQR4SFW3BSVKYwqBuDjYs2c3iYmRhISk0KdPT2rV2kKZMo2xnGEKA79jHKqaBgxU1e2qOlxVh1nSMCZrZ86cYeDA/1G//ixCQlIACAlJoWHDWWzc6G1sxuSWQAbHF4rIw2Lta2P8Wrt2LTVq1GDp0kHMm9eY5GSnQZ+cXJz4+ObY1bamsAgkcbwITAV+szvHjfmrU6dO0bVrV2rXrs2RI0eYObMtN9xwI+HhX9Kt21tERiZRrlxZGjXyOlJjcoffwXF3jKOOqn6edyHlPhscN8GyfPlyoqOj2bVrF9WrV2f9+vXunvSrqiKAjfheVdWrVy9iYmLyPlhjcuisBsdVNU1E3sG5AdAY4zp27Bjdu3dn5MiRVKhQgcWLF3PPPfd4HZYxecLGOIzJobi4OMLDwxk1ahSdOnUiMTHRkoYpUgK5AfBFoASQIiKncS7JVVW9JKiRGZPP/Pzzz3Tq1IkPP/yQSpUqsWrVKmrXru11WMbkOZsd15gATJs2jeeff56ff/6ZV199lVdffZXzzjvP67CM8US2iUNE7sxsu6ouz/1wjMlfDh48yPPPP8+MGTOoXr06CxcupGrVql6HZYynAumq6urz/HygFrAOsE5dU2ilT33+n//8h1OnTtG3b1+6dOlC8eKB/C9jTOEWSFdVE9/XInIN8HbQIjLGY7t376Zt27YsWLCAunXrEhsby8033+x1WMbkG4FcVZXRXiA8twMxxmtpaWmMGDGCsLAwVqxYwX//+1+WL19uScOYDAIZ4/gv7pTqOIkmAtgUxJiMyXM7d+6kdevWrFixgvr16zNmzBjKly/vdVjG5EuBtDgScMY01gFfAC+p6hOBHFxEGorIDhHZJSLdM9kvIjLU3b9ZRKrnoG4XEVERuSKQWIzJTEpKCgMGDKBq1aokJiYyfvx44uPjLWkY40eWLQ53wabSqvpBhu1hIlJaVQ/7O7C7+NNwoD5O99ZaEZmjqlt9ijUCKrqPSGAkEJldXXecpT6wO0ef1hgfmzdvJjo6moSEBB566CFGjBhBmTJlvA7LmHzPX4vjv0DpTLaXA4YEcOxawC5V/UZVzwCTgWYZyjQDJqhjNXCpiJQJoO67OGuE+F+FyphMnDlzhl69elGjRg12797NlClTmDFjhiUNYwLkL3FUVtVlGTeqajxQJYBjlwX2+Lze624LpEyWdUWkKbBPVf2Os4jIMyKSICIJhw/7bRyZIuTLL7+kevXqvPHGG7Ro0YKtW7fyz3/+01blMyYH/CWOkLPcly6z/xMzthCyKpPpdhG5EOgB9MzuzVV1jKrWVNWapUtn1nAyRcmpU6fo0qULderU4ejRo8ybN48PP/yQUqVKeR2aMQWOv8TxlYg0zrhRRBoB3wRw7L3ANT6vywH7AyyT1fYbgArAJhH5zt2+XkSuDiAeU0QtW7aMKlWqMHDgQNq0aUNSUhKNG//l1DbGBMjf5bj/AeaKyL9wrqgCqIkzxfqDARx7LVBRRCoA+4AWwGMZyswB2ovIZJzB8aOqekBEDmdWV1WTgCvTK7vJo6aq/hhAPKaIOXbsGC+99BKjRo3ihhtuYMmSJdSrV8/rsIwp8LJMHKq6U0Qq4/yxT7/hbxnQVlVPZ3dgVU0RkfZAPFAMGK+qSSLyrLt/FDAfaAzsAk4B//ZX9yw/oymC5s+fT9u2bdm/fz+dO3fmjTfe4MILL/Q6LGMKBb8rABYWtgJg0fHTTz/RqVMnJk6cSFhYGLGxsURGRnodljEFUlYrAJ7NlCPG5DuqytSpUwkNDWXy5Mn07NmTdevWWdIwJghsqk9T4B04cIDnn3+emTNnUqNGDRYtWkSVKoFcMW6MORsBtzhEJEREqonIldmXNib4VJX33nuP0NBQ4uLi6N+/P6tXr7akYUyQZZk4RGSUiIS5z0viTGw4AdggIi3zKD5jMvX999/TsGFDoqKiqFy5Mps2baJbt262XoYxecBfi+MOnyuZ/g3sVNXKQA2c6T6MyXNpaWkMGzaMsLAwVq1axfDhw1m6dCk33XST16EZU2T4+3l2xud5fWAqgKoetOkZjBd27NhB69atWblyJffffz+jR4/muuuu8zosY4ocfy2OX0TkQRGpBtQFFgCISHHggrwIzhhwpj7v378/VatWJSkpiffff5+4uDhLGsZ4xF+Loy0wFLga6KSqB93t9wLzgh2YMQCbNm0iKiqK9evX83//938MHz6cq6+2GWaM8ZLfO8eBhplsj8e5o9uYoPntt9/o3bs3/fr14/LLL2fq1Kk88sgjXodljMH/Qk5D/VVU1Y65H44xsHr1aqKjo9m6dStPPfUUgwYNsllsjclH/HVVPQtsAabgzExrI+ImqE6ePMlrr73G4MGDKVu2LPPnz6dRo0Zeh2WMycBf4igD/BN4FEgBPgGmq+qRvAjMFC1LliyhdevWfPPNN7Rr145+/fpxySWXeB2WMSYTWV5Vpao/qeooVb0beBq4FEgSkSfzKDZTBBw9epRnnnmGe+65h3POOYelS5cyYsQISxrG5GPZ3mYrItWBljj3csTxx9ocxuRIairExcGGDVCtGqSmzuX555/lwIEDdOnShddff92mPjemAPA35cjrIrIOeBFnHY6aqhqtqlvzLDpT4MXExCAiiBTjkkvm0LnzOo4f70vnzut47LE09u07QFpaGiVKlLCkYUwBkeV6HCKShrNE7K/upvSCAqiqFpiZ5Gw9Du/NnQu9eu1g9epwQkJSSE4uTrVqG3jzzZtp3jyQJeyNMXktq/U4/HVVVQhiPKaIWbbsGPfdN5OQkBQAQkJSaNJkHklJ4TRv7nFwxpgc8Tc4/n1mD2AvcHvehWgKMlUlNjaWESPaMGtWfZKTnd8qycnFiY9vTkSEt/EZY3LO3xjHJSLysogME5EG4uiA0331r7wL0RRU3377LQ0aNKB169bUqHGIa665lvDwL+nW7S0iI5MoV64sdpuGMQWPvzGO2cAR4Auc+akuA84FXlDVjXkVYG6wMY68lT71+csvv0xKSgpnzqRPtHwO0AiIADbiXKSXBkCvXr2IiYnJ+2CNMVk6mzGO6931NxCRccCPwLWqejxIMZpCYPv27URHR7Nq1SoaNmzI6NGjufbaa70OyxiTi/xNq56c/kRVU4FvLWmYrCQnJ9O3b18iIiLYtm0bH3zwAfPnz7ekYUwh5K/FUVVEjrnPBbjAfZ1+Oa7d2msA2LhxI1FRUWzYsIGHH36YYcOG2dTnxhRi/qZVL5aXgZiC5/Tp0/Tu3Zv+/ftTqlQppk2bxsMPP+x1WMaYIMt2yhFjMvPFF18QFRXF9u3badWqFYMGDeLyyy/3OixjTB7wN8ZhzF+cPHmSTp06UbduXU6dOkVcXBzvv/++JQ1jihBrcZiALV68mDZt2vDtt9/y/PPP07dvXy6++GKvwzLG5DFrcZhsHT16lDZt2nDfffdRvHhxli1bxrBhwyxpGFNEWeIwfn366aeEhoYyfvx4unXrxqZNm7jzzju9DssY4yFLHCZThw8f5rHHHqNp06aUKlWKNWvW0L9/fy644AKvQzPGeMwSh/kTVWXy5MmEhoYybdo0Xn/9dRISEqhZ8y+zDhhjiigbHDe/279/P+3atWPOnDnUqlWL8ePHExYW5nVYxph8xloc5vepz0NDQ1m0aBHvvPMOq1atsqRhjMmUtTiKuG+//ZY2bdqwePFi6tWrx9ixY7nxxhu9DssYk48FtcUhIg1FZIeI7BKR7pnsFxEZ6u7fLCLVs6srIgNEZLtbfqaIXBrMz1BYpaamMmTIEMLDw/nyyy8ZNWoUixcvtqRhjMlW0BKHiBQDhuMswBAKtBSR0AzFGgEV3cczwMgA6i4Cwt01z3cCLwfrMxRW27Zt484776RTp07Uq1ePpKQk2rZtyznnWM+lMSZ7wfxLUQvYparfqOoZYDLQLEOZZsAEdawGLhWRMv7qqupCVU1x668GygXxMxRoqakwdy68+abz7+nTybz11ltERESwfft2Jk6cyNy5c7nmmmu8DtUYU4AEM3GUBfb4vN7rbgukTCB1AaJwlpH7CxF5RkQSRCTh8OHDOQy94IqJiUFEECnGJZfMoXPndRw/3pfOnddRqlQcPXq8xpkzZ4iKiuLxxx9HRLwO2RhTwAQzcWT2FynjOrVZlcm2roj0AFKAjzJ7c1Udo6o1VbVm6dKlAwi3cIiJiUFV+fTTVG655Wa2bKnN22+/wpYttSlX7hpeeWUlqsqAAQO8DtUYU0AFM3HsBXz7QMoB+wMs47euiLQCHgQe16wWTS/iNmyA+vVnERLi9OqFhKTQvPlCLrigjseRGWMKumAmjrVARRGpICLnAi2AORnKzAGecq+uqg0cVdUD/uqKSEPgJaCpqp4KYvwF1smTJ1m7dgwzZ95HcrJzxXVycnEWLmxORIS3sRljCr6gJQ53ALs9EA9sA6aoapKIPCsiz7rF5gPfALuAscBz/uq6dYYBFwOLRGSjiIwK1mcoiJ566ikuuugiPv20HXv37iMsbA3dur1FWNgaduzYTpMmxRARYmJivA7VGFNASVHo6alZs6YmJCR4HUZQ/fLLL3Tp0oXY2FhuuukmYmNjqVPnduLiYONGiIiARo2gmC0IbIwJkIisU9W/TFRnd44XAnPmzKFdu3b88MMPdO/enZ49e/4+i+2DDzoPY4zJLXbHVwF2+PBhWrZsSbNmzbjiiitYs2YNffv2tanPjTFBZYmjAFJVJk2aRGhoKDNmzODNN99k7dq11KhRw+vQjDFFgHVVFTD79u3j2WefZe7cuURGRjJ+/HhCQzPO5GKMMcFjLY4CQlUZO3YsoaGhLF68mEGDBvH5559b0jDG5DlrcRQAX3/9NW3atGHJkiXcfffdjB07lhtuuMHrsIwxRZS1OPKx1NRU3n33XSpXrkxCQgKjR49m8eLFljSMMZ6yFkc+tXXrVqKjo1m9ejUPPPAAo0aNolw5mwjYGOM9a3HkM8nJyfTu3Ztq1arx1VdfMXHiRD799FNLGsaYfMNaHPnI+vXriYqKYtOmTfzrX//iv//9L1deeaXXYRljzJ9YiyMfOH36NC+//DK1atXi0KFDzJw5k08++cSShjEmX7IWh8dWrlxJdHQ0O3fuJCoqinfeeYfLLrvM67CMMSZL1uLwyIkTJ+jQoQN33nknv/32GwsXLiQ2NtaShjEm37PE4YFFixYRHh7O8OHDad++PVu2bKF+/fpeh2WMMQGxxJGHjhw5QlRUFA0aNOD8889nxYoVDB06lIsuusjr0IwxJmCWOPLIrFmzCA0NZcKECXTv3p2NGzdSt25dr8Myxpgcs8HxIDt06BAdOnRgypQpVK1alXnz5lG9enWvwzLGmLNmLY4gUVU++ugjQkNDmTVrFr1792bt2rWWNIwxBZ61OIJgz549tGvXjnnz5lG7dm1iY2NtFltjTKFhLY5clJaWxujRowkLC2PJkiUMGjSIlStXWtIwxhQq1uLIJV9//TWtW7dm6dKl3HPPPYwdO5brr7/e67CMMSbXWYvjb0pNTWXQoEFUrlyZ9evXM2bMGD777DNLGsaYQstaHH9DUlIS0dHRrFmzhgcffJCRI0faLLbGmELPWhxn4cyZM7z55ptUq1aNXbt28fHHHzNnzhxLGsaYIsFaHDmUkJBAdHQ0mzdvpkWLFgwZMsRmsTXGFCmWOPxITYW4ONiwAUJDf2P16l4MGjSAq666itmzZ9O0aVOvQzTGmDxniSMLqanQvPlJ9u3bS/36s+jRowF79txGq1ZRDBo0gEsvvdTrEI0xxhM2xpFBTEwMIkLx4g+yY8d2Vq8Op1+/7iQm1qJs2XK8994BLrvsMmJiYrwO1RhjPGEtjgxiYmKIiYmhVatdXHXVVEJCUgAICUnh4YcXUaLEXF591eMgjTHGQ9biyOCll15CRJgwoROzZzcgOdnJrcnJxZk+vT6vvfYgImItDmNMkWUtDh8zZsxgwoQJFCtWjK5dI0hMvInIyCTuv38m8fHNueWWsmzbNpdixbyO1BhjvGOJA/jhhx9o374906ZNIyIigvnz51OtWjX3qqqL2bjxJd54Axo1wpKGMabIK9KJQ1WZOHEinTp14sSJE/Tp04euXbsSEhICOEniwQedhzHGGEeRTRx79uyhbdu2xMXFUadOHWJjY6lUqZLXYRljTL5XZBLH3LnOjXxVq6axd+9YunfvSmpqKoMHD6Z9+/YUsz4oY4wJSFCvqhKRhiKyQ0R2iUj3TPaLiAx1928WkerZ1RWRy0VkkYh85f57WXZx7NqVRq9eOzh5sj/du2+ma9cy3HprbbZs2cILL7xgScMYY3JAVDU4BxYpBuwE6gN7gbVAS1Xd6lOmMdABaAxEAkNUNdJfXRF5G/hZVfu5CeUyVX3JXywlSoTrL7/sICQkheTk4lSvvpG33gqlSRPJ/Q9ujDGFhIisU9WaGbcHs8VRC9ilqt+o6hlgMtAsQ5lmwAR1rAYuFZEy2dRtBnzgPv8AeCi7QC6++Jc/3cjXqNEcmjZ9DRGhXr16f+tDGmNMURPMMY6ywB6f13txWhXZlSmbTd2rVPUAgKoeEJFMp6YVkWeAZwDOOackNWqACKhCUtKHaXDwG+DosmXLEPGs5XEF8KNXb+6HxZUzFlfOWFw542Vc12W2MZiJI7O/xhn7xbIqE0hdv1R1DDAGQEQS1q37a3PLayKSkFkz0GsWV85YXDljceVMfowrmF1Ve4FrfF6XA/YHWMZf3R/c7izcfw/lYszGGGOyEczEsRaoKCIVRORcoAUwJ0OZOcBT7tVVtYGjbjeUv7pzgFbu81bA7CB+BmOMMRkEratKVVNEpD0QDxQDxqtqkog86+4fBczHuaJqF3AK+Le/uu6h+wFTRCQa2A38M4BwxuTeJ8tVFlfOWFw5Y3HljMUVoKBdjmuMMaZwsmnVjTHG5IglDmOMMTlSoBNHMKY0yaO4Hnfj2Swiq0Skqs++70QkUUQ2ikhCHsdVT0SOuu+9UUR6Blo3D2Lr6hPXFhFJFZHL3X1B+c5EZLyIHBKRLVns9+r8yi4ur86v7OLy5PwKIK48P7fcY18jIktEZJuIJInIC5mU8eQcy5aqFsgHzqD518D1wLnAJiA0Q5nGQBzOfSG1gTWB1g1yXLfhTJUC0Cg9Lvf1d8AVHn1f9YC5Z1M32LFlKN8E+F8efGd3AtWBLVnsz/PzK8C48vz8CjAur84vv3F5cW65xy4DVHefX4wzzZLnf8MCeRTkFkewpjQJelyqukpVj7gvV+PcpxJsf+czB/P7OpvjtwQm5eL7Z0pVlwM/+ynixfmVbVwenV+BfF9Z8fT7yiBPzi1wZr5Q1fXu8+PANpxZM3x5co5lpyAnjqymKwmkTCB1gxmXr2icXxTpFFgoIuvEmTYltwQaVx0R2SQicSISlsO6wY4NEbkQaAhM99kcrO8sO16cXzmVV+dXoLw4vwLi5bklIuWBasCaDLvy5TlWkNfj8HRKEz8CPraI3I3zP/btPpvrqup+cebgWiQi291fTHkR13rgOlU9Ic7MxbOAigHWDXZs6ZoAn6uq7y/IYH1n2fHi/ApYHp9fgfDq/AqUJ+eWiFyEk6w6qeqxjLszqeL5OVaQWxzBmtIkL+JCRKoA44BmqvpT+nZV3e/+ewiYidMkzZO4VPWYqp5wn88HQkTkikDqBjs2Hy3I0JUQxO8sO16cXwHx4PzKlofnV6Dy/NwSkRCcpPGRqs7IpEj+PMfyajAltx84raVvgAr8MTgUlqHMA/x5YOnLQOsGOa5rce6Wvy3D9hLAxT7PVwEN8zCuq/njptBaOHfmSzC/r5z89wBK4vRVl8iL78w9ZnmyHuzN8/MrwLjy/PwKMC5Pzq/s4vLw3BJgAjDYTxnPzjF/jwLbVaXBm9IkL+LqCZQCRogzpXuKOrNfXgXMdLcVBz5W1QV5GNcjQDsRSQF+BVqoc5YG7fvKQWwAzYGFqnrSp3rQvjMRmYRzJdAVIrIX6AWE+MSU5+dXgHHl+fkVYFyenF8BxAV5fG656gJPAokistHd9gpO4vf0HMuOTTlijDEmRwryGIcxxhgPWOIwxhiTI5Y4jDHG5IglDmOMMTliicMYY0yOWOIwBZ47m2n6zKafisil7vbyIqIi8qZP2StEJFlEhrmvbxaRpW79bSLyl9XWRORbEbk5w7bBItLNT0wnfGJ4LJc+ql8i0klEnjrLuktFpKb7/JUAyk8WkYpn816m4LPEYQqDX1U1QlXDcW7iet5n3zfAgz6v/wn4Xu8+FHjXrV8J+G8mx5+Mc1cxACJyDs49CZ8EEFt5IOiJQ0SKA1HAx7lwuGwTBzASyDJxmsLNEocpbL7gz5O9/QpsS/81DTwKTPHZXwZn+gYAVDUxk2NOwidx4EzT/Z2qfi8iL7otnS0i0imTuv2AO9wWzX/cFsgKEVnvPm4DJxmJyAhx1mWYKyLzReQRd18NEVnmTrQX786OmtE9wHr3xrBKIvJl+g73PTe7z+8VkQ3irDExXkTO8z2IiPQDLnDj/UhESojIPHdiwi0i8qhbdAVwn5uwTBFjicMUGiJSDLgXmJNh12SghYiUA1L585w+7wL/E2e21v+kd3P5UtXNQJr8sSBSC2CSiNTAuZM3Emc6iDYiUi1D9e7ACrdF8y5wCKivqtVxkthQt9z/4bROKgOtgTruZwrBaQU9oqo1gPFAn0w+fl1gnRvvNuBcEbne3fcoMEVEzgfeBx5V1co4d0O3y/BZu/NHC+5xnNli96tqVbdFt8Atl4ZzN3NVTJFjicMUBhe4Uzb8BFwOLMqwfwFQH2ethT91L6nqe0AlYCrOtBSrM/4Kd03CST7FcdY9mIoz6+xMVT2pzuR9M4A7sok1BBgrIonuMULd7bcDU1U1TVUPAkvc7TcD4Tgzs24EXiXz9TXKAId9Xk8B/uU+f9T93DcD36rqTnf7BzitJ38ScVoW/UXkDlU96rPvEPCPbOqbQsgShykMflXVCOA6nAnffMc4UGehm3VAZ/681kL6/v2qOl5VmwEpOH+oM5qE84f4PmCzOrOlZja1dXb+A/yA80u9phsvfo4lQJLbAohQ1cqq2iCTcr8C5/u8/gT4l4jcBKiqfnU28bpJpgZOAukrPsu9uu/3a06PaQo+Sxym0HB/DXcEurhdPL4GAi+pzxTj8Pu6zSHu86txJgfcl8mxv8Zp0fTjj6m3lwMPiciFIlICZ6K8FRmqHsdZFjRdSeCA29XzJM4EdQArgYfdsY6rcFo/ADuA0iLye9eV/LEAkq9twI0Z4k0FXuOPVtZ2oLyIpJd7EliWybGSfb6TfwCnVHUi8A7OEqzpbuLPFxqYIsIGtkyhoqobRGQTzjjECp/tSWT+R64BMERETruvu7pdRZmZBPTFWZcBVV0vIu8D6QPR41R1Q4Y6m3Fmf92EM74wApguIv/E6Y5Kn411Os74zBactafXAEdV9Yw7SD5UREri/D87OJPPEgd8mGHbJ8AAnKm3UdXTIvJvYKrb5bYWGMVfjQE2i8h6nGm/B4hIGpCMOybiJrdfVfVAFt+VKcRsdlxj8gkRuUid1fFK4SSjun6SWGb1ZwLd3G6poBKR/wDHVDU22O9l8h9rcRiTf8x1r+o6F3gzJ0nD1R1nkDzoiQP4hb+2cEwRYS0OY4wxOWKD48YYY3LEEocxxpgcscRhjDEmRyxxGGOMyRFLHMYYY3Lk/wFB2xUqMTvOiQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Create a range of $x$ values for plotting:\n", "xx = np.arange(-1, 10, 0.1)\n", "\n", "# Plot of the data as we did above...\n", "plt.errorbar(Vrms, Irms, errIrms, fmt = 'ko', markersize = 5,\\\n", " linewidth = 1.8,\\\n", " markeredgecolor = 'b',\\\n", " markerfacecolor = 'yellow',\\\n", " capsize = 5)\n", "\n", "# Plot the best-fit line...\n", "plt.plot(xx, fitFcn(xx), 'k-')\n", "\n", "# Add some plot labels...\n", "plt.xlabel('RMS Voltage (volts)')\n", "plt.ylabel('RMS Current (amps)')\n", "plt.title('Unweighted Linear Fit')\n", "\n", "# Set the range of the axes...\n", "plt.axis((0, 2.2, 0, 0.01));" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "All of this has produced an \"unweighted\" fit to the data. To include weights, all we need to do is include another option in *curve_fit()*. Everything else is exactly the same! The new option is *sigma* and it is simply a list of the errors in the $y$-values. Note that many fitting routines require you to provide the actual weights as $1/\\sigma^2$. That is not the case here. You just have to provide the absolute $y$-uncertainties." ] }, { "cell_type": "code", "execution_count": 85, "metadata": {}, "outputs": [], "source": [ "a_fit, cov = curve_fit(linearFunc, Vrms, Irms, sigma = errIrms)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Extract the fit parameters and their uncertainties." ] }, { "cell_type": "code", "execution_count": 86, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The best-fit parameters are: (1) slope = 0.004560548503316531 and (2) intercept -2.8519759527912626e-05\n", "The error in the slope is: 9.344743968455787e-06 \n", "The error in the intercept is: 1.7648857220463362e-06\n" ] } ], "source": [ "print('The best-fit parameters are: (1) slope =', a_fit[1], 'and (2) intercept',\\\n", " a_fit[0])\n", "print('The error in the slope is:', np.sqrt(np.diag(cov))[1], '\\n'\n", " 'The error in the intercept is:', np.sqrt(np.diag(cov))[0])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Get the best-fit line." ] }, { "cell_type": "code", "execution_count": 87, "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$x \\mapsto \\text{-2.8519759527912626e-05} + \\text{0.004560548503316531}\\,x$" ], "text/plain": [ "Polynomial([-2.85197595e-05, 4.56054850e-03], domain=[-1, 1], window=[-1, 1])" ] }, "execution_count": 87, "metadata": {}, "output_type": "execute_result" } ], "source": [ "fitFcn = np.polynomial.Polynomial(a_fit)\n", "fitFcn" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Plot the data..." ] }, { "cell_type": "code", "execution_count": 88, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEWCAYAAABxMXBSAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA32ElEQVR4nO3dd3gU5drH8e8tBAtYeBUVAcWC4hJ6JIpipRjEgxy7ciwJ9QCKIAiiElRUmgpIEQiKiLRDlY5IVRAChBKQIiIgIKL0mnK/f8yga0w2G8xmUu7Pde3F7rT97Trm3meemWdEVTHGGGOCdY7XAYwxxuQtVjiMMcZkiRUOY4wxWWKFwxhjTJZY4TDGGJMlVjiMMcZkiRUOk++JyNMiMifIZZ8TkSUhzHJW2xeRV0VkWCgyZRcRSRSRu73OYULPCofJlUSks4jMSDNtSwbTngi0LVUdpap1synXAhFpkh3bSmfbZUVERaRw2nmq+o6qhuR9s8r9Dk6KyFG/x22qWkFVF7jLxIrI5x5HNSFihcPkVouA20WkEICIXAmEAdXSTLvBXdaEwJnvOh2tVbWY32NpjgYznrLCYXKrFTiFoor7+k5gPrApzbQfVHW3iFwsInEiskdEfhaRt/0KzF8OD4lIXRHZJCKHRGSgiCxM24oQkd4ickBEfhSRKHdad6AW8JH7K/sjd3p5EZkrIr+7233MbzuXishUETksIsuB68/my/D/Be/XMnlWRHaIyH4R6eK37Dki0klEfhCR30RknIj8n9/88SKy1/38i0Skgt+8T0VkkIjMEJFjwD1ZyLhdRGqLyP3Aq8Dj7ve05mw+s8m9rHCYXElVTwPf4RQH3H8XA0vSTDvT2hgBJOO0QKoCdYG/HdoRkcuA/wGdgUtxClHNNItFutMvA3oCcSIiqtrFzXDm13ZrESkKzAW+AC4HngQG+v0xHgCcBEoC0e4ju9wB3ATcB7whIje7018AHgLuAq4CDrg5zpgJlHPzrgJGpdnuU0B34EKc7ztLVHUW8A4w1v2eKmd1GyZ3s8JhcrOF/FkkauH80V6cZtpCEbkCiALaquoxVd0HfACk1/dRH0hU1Ymqmgz0A/amWeYnVR2qqik4BakkcEUGGRsA21X1E1VNVtVVwATgEbfF8zDwhptrvbu97NJNVU+o6hpgDXDmD3RzoIuq7lLVU0Csm6cwgKoOV9UjfvMqi8jFftudoqrfqGqqqp7M4L37ichB97EqGz+TyQP+1glnTC6yCGglIsWBEqq6RUR+AUa408LdZa7BOay1R0TOrHsOsDOdbV7lP11VVUR2pVlmr9/84+42i2WQ8RogUkQO+k0rDIwESrjP/XP8lOGnzTr/gnfcL+M1wCQRSfWbnwJcISJ7cVoTj7r5zixzGXDIfZ7e95bWC6qaq8/yMqFjhcPkZkuBi4FmwDcAqnpYRHa703ar6o8ichI4BVzmtiIC2QOUPvNCnKpQOuPF/ybtcNI7gYWqWiftgm6LIxkoA3zvTr46C+91tnYC0ar6TTqZ/gM0BGoD23G+3wOA+C2WHUNm27Db+ZgdqjK5lqqeAOKBdjiHqM5Y4k5b5C63B5gD9BGRi9zO4etF5K50NjsdqCgiD7mHbloBV2Yh1i/AdX6vpwE3ish/RCTMfdwiIje7h7omArEicoGI+IBng3iPc0XkPL9HVv8/HQx0F5FrAESkhIg0dOddiFNkfwMuwOmLCIVfgLJnkd3kAfYf1eR2C3E6cf07aRe70/xPw30GKAJswPkF/T+cvom/UNX9OIdpeuL88fThFKdTQebpi9NfcEBE+qnqEZyO+CeA3TiHj3oA57rLt8Y5hLQX+BT4JIj3OAqc8HvcG2Q2/4xTgTkicgRYhtPhD/AZzuGyn3G+q2VZ3Hawxrv//mZ9IPmP2I2cTEHm/iLeBTytqvO9zmNMXmAtDlPgiEg9EblERM7Fud5ACN0vb2PynZAWDhG5370gaquIdEpnvohIP3f+WhGp5jdvuIjsE5H1adb5P/diqy3uv8VD+RlMvnQb8AOwH3gQeMjtTzHGBCFkh6rcM0o2A3VwDgWsAJ5U1Q1+y9QH2uCcWx8J9FXVSHfenTjHej9T1XC/dXoCv6vqe24xKq6qr4TkQxhjjPmbULY4agBbVXWbexXwGJzTAP01xCkMqqrLgEtEpCSAqi4Cfk9nuw358yKqEThXyBpjjMkhobyOoxR/vZBoF3+e2RFomVI459pn5Ar39EtUdY+IXJ7eQiLSDOdcf4oWLVq9fPnyWUtvjDEF0PHjx9m+fTsnTpwA2K+qJdIuE8rCIelMS3tcLJhlzoqqDgGGAERERGh8fHx2bNYYY/KlkydP8uabb9KzZ09KlCjBqFGj+Pe//53uSAehLBy7cK6YPaM0znnuWV0mrV9EpKTb2igJ7PvHSY0xpgD79ttviYmJ4fvvv+e5557j/fffp3jxjM87CmUfxwqgnIhcKyJFcC6QmppmmanAM+7ZVbcCh84chgpgKn9effssMCU7QxtjTEFx9OhRXnzxRe644w6OHz/OrFmz+OSTTwIWDQhh4XDHDGoNzAY2AuNUNVFEWohIC3exGcA2YCswFPjvmfVFZDTOWEU3icguEYlxZ70H1BGRLThnbL0Xqs9gjDH51dy5c6lYsSL9+vWjVatWrF+/nnr16gW1boG4ctz6OIwxxnHw4EHat2/P8OHDufHGG4mLi+OOO+5Id1kRWamqEWmn25XjxhhTQEyePBmfz8eIESPo1KkTCQkJGRaNQGxYdWOMyef27dtHmzZtGDduHJUrV+bLL7+kevXqZ709a3EYY0w+paqMGjUKn8/H5MmTeeutt1ixYsU/KhpgLQ5jjMmXdu3aRYsWLZg+fTqRkZEMHz4cn8+XLdu2FocxxuQjqampfPzxx/h8Pr7++mvef/99vvnmm2wrGmAtDmOMyTe2bt1K06ZNWbBgAffeey9Dhw7luuuuy3zFLLIWhzHG5HEpKSn06dOHSpUqsWrVKoYMGcJXX30VkqIB1uIwxpg8LTExkejoaJYvX06DBg0YNGgQpUuXDul7WovDGGPyoNOnT/Pmm29StWpVtm3bxhdffMHUqVNDXjTAWhzGGJPnxMfHEx0dzbp163jiiSfo168fJUr8bfTzkLEWhzHGeCw2NhYRyfTRpUsXOnbsSGRkJL/99htTpkxh9OjROVo0wMaqMsaYXOeqq64CYPfuP+8ysWjRIpo0acKWLVto0qQJvXr14pJLLglpDhuryhhj8oCUFDh5sjZHjrRl2jQ4ePAIrVq14q677iI5OZmvvvqKoUOHhrxoBGJ9HMYYk0ukpECjese4/GBb/qWzee3RDezQXRw49TFt27bl7bffpmjRol7HtMJhjDG5xcyZ8PN3u1inkYSRTPeTb1BRVvBqr3W8/PLNXsf7gx2qMsaYXGL1aqhzbDJhJAMQRjIPMZuTJ3NP0QArHMYYkyvs3buXmTPfYZLWJsk9GJREYeYUbUSVKt5mS8sKhzHGeEhVGTFiBD6fj5Ur3+bc688nnOV05B0iiyVSOrIUUVFep/wrKxzGGOORHTt2UL9+fZ577jl8Ph9r165m9SYfvxb/gEHFzuHN0TcyaXZRChXyOulfWeEwxpgclpqaysCBA6lQoQKLFy8mKiqKb775hvLly1O4sHDgwEiOHu3Egw8KhQv/eQFgbGys19EBO6vKGGNy1ObNm2nSpAmLFy+mTp06DBkyhLJly3odK0usxWGMMTkgOTmZnj17UrlyZdatW8fw4cOZPXt2nisaYC0OY4wJubVr1xITE0N8fDwPPfQQAwcOpGTJkl7HOmvW4jDGmBA5deoUXbt2pXr16uzYsYNx48YxceLEPF00wFocxhgTEt999x0xMTEkJibSuHFjPvzwQy699FKvY2ULa3EYY0w2On78OO3bt6dmzZocOnSIadOmMXLkyHxTNMBaHMYYk20WLFhAkyZN+OGHH2jRogU9evTgoosu8jpWtrMWhzHG/EOHDh2iefPm3HPPPQDMnz+fQYMG5cuiAVY4jDHmH5k+fToVKlRg2LBhtG/fnrVr13L33Xd7HSukrHAYY8xZ2L9/P40bN6ZBgwZccsklLF26lN69e3PBBRd4HS3krHAYY0wWqCrjxo3D5/MxduxYunbtyqpVq6hRo4bX0XKMdY4bY0yQ9uzZw3//+18mT55MREQE8+bNo2LFil7HynHW4jDGmEyoKsOHD+fmm29m1qxZ9OzZk6VLlxbIogHW4jDGmIC2b99Os2bNmDt3LrVq1WLYsGHceOONXsfylLU4jDEmHampqfTr14/w8HCWLl3KwIEDWbBgQYEvGhDiwiEi94vIJhHZKiKd0pkvItLPnb9WRKpltq6IVBGRZSKSICLxIlJweqSMMdkiNjb2j3tcZPQoVKgQL774IrVq1SIxMZGWLVtyzjn2WxtCWDhEpBAwAIgCfMCTIuJLs1gUUM59NAMGBbFuT6CbqlYB3nBfG2NM0GJjY1FVVJXkZKV48f9QrFgPJk1K4u233+Xcc8+lePHijBgxghkzZnD11Vd7HTlXCWUfRw1gq6puAxCRMUBDYIPfMg2Bz1RVgWUicomIlATKBlhXgTOXY14M7A7hZzDG5GMpKdCo3jEuP9iWf+lsOj+8jh2pPh74978YMKAfV155pdcRc6VQFo5SwE6/17uAyCCWKZXJum2B2SLSG6fFVDO9NxeRZjitGPu1YIxJ18yZ8POyXazTSMJIpnvqG1Q/bw3PPz8OqxkZC+UBO0lnmga5TKB1WwIvqWoZ4CUgLr03V9UhqhqhqhElSpQIMrIxpiCZOnUH9x2bSBjJAISRzAOnviQhwdtcuV0oC8cuoIzf69L8/bBSRssEWvdZYKL7fDzOITFjjAnasWPHaNu2LUOHtmKq1CPJPfiSRGFmF21ElSre5svtQlk4VgDlRORaESkCPAFMTbPMVOAZ9+yqW4FDqronk3V3A3e5z+8FtoTwMxhj8pmvvvqK8PBw+vbtS8uWZbnuzhsIl+V05B0iiyVSOrIUUVFep8zdQtbHoarJItIamA0UAoaraqKItHDnDwZmAPWBrcBx4PlA67qbbgr0FZHCwEncfgxjjAnk4MGDvPzyy8TFxVGuXDkWLVpErVq1SEmBEiU+YFBSBUaPvpGoKChUyOu0uZs4JzTlbxERERofH+91DGOMR6ZOnUrLli3Zu3cvL7/8MoULF+add97JdL2uXbsSGxsb+oC5lIisVNWItNNtyBFjTL7166+/8sILLzBmzBgqVqzIlClTiIhw/g52797d43R5l10GaYzJd1SV0aNH4/P5mDBhAt26dSM+Pv6PomH+GWtxGGPylZ9//pkWLVowbdo0atSowfDhw6lQoYLXsfIVa3EYY/IFVWXo0KH4fD7mzZtHnz59+Pbbb61ohIC1OIwxed4PP/xA06ZNmT9/PnfffTdDhw7lhhtu8DpWvmUtDmNMnpWSksIHH3xAxYoViY+PZ/DgwcybN8+KRohZi8MYkydt2LCBmJgYli1bRv369Rk8eDBlypTJfEXzj1mLwxiTpyQlJfH2229TtWpVtmzZwueff860adOsaOQga3EYY/KMVatWER0dzZo1a3jsscfo378/l19+udexChxrcRhjcr0TJ07QuXNnatSowb59+5g0aRJjx461ouERa3EYY3K1JUuWEBMTw+bNm3n++efp06cPxYsX9zpWgRawcIjIeUADoBZwFXACWA9M9xt00Bhjst2RI0d49dVXGTBgAFdffTVz5syhTp06XscyBCgcIhILPAgsAL4D9gHnATcC77lFpb2qrg19TGNMQTJnzhyaNm3Kzp07ad26Ne+88w7FihXzOpZxBWpxrFDV2AzmvS8ilwN2T1ZjTLY5cOAA7dq149NPP+Wmm25i8eLF3H777V7HMmlk2DmuqtPTThORc0TkInf+PlW1scqNMdli4sSJ+Hw+Ro4cSefOnUlISLCikUtlelaViHwhIheJSFFgA7BJRDqEPpoxpiD45ZdfePTRR3n44Ye58sorWb58Oe+88w7nnXee19FMBoI5HdenqoeBh3Du2Hc18J9QhjLG5H+qysiRI/H5fEydOpXu3buzfPlyqlWr5nU0k4lgCkeYiIThFI4pqpoE5P/bBhpjQmbHjh088MADPPPMM9x0000kJCTw6quvEhYW5nU0E4RgCsfHwHagKLBIRK4BDocylDEmf0pNTWXQoEFUqFCBhQsX8uGHH7J48WJuvvlmr6OZLMj0AkBV7Qf085v0k4jcE7pIxpj8aMuWLTRt2pSFCxdy3333MXToUK699lqvY5mzEEzn+KUi0k9EVonIShHpC1ycA9mMMflAcnIyvXv3plKlSiQkJDBs2DDmzp1rRSMPC+ZQ1RjgV+Bh4BH3+dhQhjLG5A/r1q2jZs2adOjQgbp16/4xFLqIeB3N/APBFI7/U9W3VPVH9/E2cEmIcxlj8rDTp08TGxtL9erV2b59O2PGjGHy5MlcddVVXkcz2SCYQQ7ni8gTwDj39SPA3y4ONMYYgBUrVhAdHc369et56qmn6Nu3L5dddpnXsUw2CqbF0Rz4AjjtPsYA7UTkiIjY2VXGGACOHz/Oyy+/zK233sqBAwf48ssvGTVqlBWNfCiYs6ouzIkgxpjcLzY2lm7durmvzgGigKrAamAmkApA48aNadCggScZTeiJaubX8olIJaAsfoVGVSeGLlb2ioiI0Ph4G1bLmOySkgKN6h1j87xN/IvZTKEuv533O6OnFqJOnXu9jmeyiYisVNWItNMzbXGIyHCgEpDImZ8TzpXjeaZwGGOy18yZsOObn1hHJGEk0503iCy0nlOnbvI6mskBwXSO36qqvpAnMcbkCb/99huvvbaAuic3E0YyAGEkU+/4ZBISXsGOUOV/wXSOLxURKxzGFHCqyvjx4/H5fKxb9xkzwxqQ5P72TKIws4s2okoVbzOanBFM4RiBUzw2ichaEVknInbXP2MKkD179vDwww/z2GOPUaZMGVaseJNr77yOcFlOR94hslgipSNLERXldVKTE4I5VDUcZxj1dfzZx2GMKQBUlREjRvDSSy9x4sQJ3nvvPdq3b0/hwoWZNBtKlPiAQUkVGD36RqKioFAhrxObnBBMi2OHqk51rxr/6cwj5MmMMZ766aefiIqK4vnnnyc8PJw1a9Zw4sQJwsLCEBEKFxYOHBjJ0aOdePBB57WI84iNjfU6vgmhTE/HFZGBOEOMfAmcOjPdTsc1Jn9KTU1l4MCBdOrUCYAePXrQsmVLzjknmN+ZJj/J6HTcYPaE83EKRl3gQfcR1HkTInK/2zeyVUQ6pTNf3JF3t7r9J9WCWVdE2rjzEkWkZzBZjDGZ27RpE3fddRdt2rTh9ttvJzExkVatWlnRMH8RzJXjz5/NhkWkEDAAqAPsAlaIyFRV3eC3WBRQzn1EAoOAyEDruvcCaQhUUtVTInL52eQzxvwpOTmZPn360LVrV84//3w++eQTnn32WRvF1qQrmAsAzwNigArAH3ePV9XoTFatAWxV1W3udsbg/MH3LxwNgc/UOV62TEQuEZGSOFepZ7RuS+A9VT3l5tgXxOc0xmRgzZo1REdHs2rVKv79738zYMAArrzySq9jmVwsmPbnSOBKoB6wECgNHAlivVLATr/Xu9xpwSwTaN0bgVoi8p2ILBSRW9J7cxFpJiLxIhL/66+/BhHXmILl1KlTvP7660RERLBr1y7Gjx/PhAkTrGiYTAVTOG5Q1deBY6o6AngAqBjEeum1cdP2xGe0TKB1CwPFgVuBDsA4Sac9rapDVDVCVSNKlCgRRFxjCo5ly5ZRrVo13n77bZ566ik2bNjAI4884nUsk0cEUziS3H8Pikg4zm1jywax3i6gjN/r0sDuIJcJtO4uYKI6luNcW2LjNhsThGPHjtGuXTtq1qzJkSNHmDFjBiNGjODSSy/1OprJQ4IpHENEpDjwGjAVp5+hRxDrrQDKici1IlIEeMJd399U4Bn37KpbgUOquieTdScD9wKIyI1AEWB/EHmMKdC+/vprKlWqxAcffECLFi1Yv349UXaptzkLwZxVNcx9ugi4LtgNq2qyiLQGZgOFgOGqmigiLdz5g4EZQH1gK3AceD7Quu6mhwPDRWQ9zo2lntVgxoY3poA6dOgQHTp0YOjQodxwww0sXLiQO++80+tYJg/L8AJAEWkMfKGq6Q4zIiLXAyVVdUkI82ULuwDQFFTTpk2jRYsW7Nmzh/bt29OtWzfOP/98r2OZPOJs7sdxKbBaRFYCK4FfcU7HvQG4C+fw0N8u6jPGeG///v28+OKLfPHFF1SsWJFJkyZxyy3pnoBoTJZlWDhUta+IfITTn3A7zs2cTgAbgf+o6o6ciWiMCZaqMnbsWNq0acOhQ4eIjY2lc+fOFClSxOtoJh8J2MehqinAXPdhjMnFdu/eTcuWLZk6dSq33HILw4cPJzw83OtYJh+yAWiMyeNUlbi4OHw+H3PmzKF3794sXbrUioYJmWDux2GMyaV+/PFHmjZtyrx587jrrrsYNmwYN9xwg9exTD6XaYtDRK4NZpoxJuekpKTQt29fwsPDWb58OYMHD+brr7+2omFyRDCHqiakM+1/2R3EGBOcjRs3UqtWLdq2bcvdd99NYmIizZs3t6HPTY7J8FCViJTHGRH3YhH5t9+si/AbJdcYkzOSkpLo1asX3bp1o1ixYnz22Wc0btzYhj43OS5QH8dNODdsugTn5k1nHAGahjCTMSaN1atXEx0dTUJCAo8++ij9+/fniiuu8DqWKaACXccxBZgiIrep6tIczGSMcZ08eZI333yTnj17UqJECSZOnEijRo28jmUKuGDOqtoqIq/ijIj7x/JB3MjJGPMPfPvtt8TExPD999/z/PPP06dPH4oXL+51LGOCKhxTgMXAV0BKaOMYY44ePUqXLl3o378/V199NbNnz6Zu3bpexzLmD8EUjgtU9ZWQJzHGMHfuXJo1a8b27dtp3bo17777LsWKFfM6ljF/Ecz5e9NEpH7IkxhTgB08eJCYmBjq1q1LkSJFWLx4Mf3797eiYXKlYArHizjF46SIHBaRIyJyONTBjCkopkyZgs/nY8SIEXTq1ImEhATuuOMOr2MZk6FgbuR0YU4EMaag2bdvHy+88AJjx46lcuXKfPnll1SvXt3rWMZkKpghR0REGovI6+7rMiJSI/TRjMmfVJVRo0bh8/mYNGkSb731FitWrLCiYfKMYA5VDQRuA55yXx8FBoQskTH52K5du3jwwQdp3Lgx5cqVY/Xq1bz22muEhYV5Hc2YoAVTOCJVtRVwEkBVDwB2VxhjsiA1NZUhQ4ZQoUIF5s+fzwcffMCSJUvw+XxeRzMmy4I5HTdJRAoBCiAiJYB070NuTEEWGxtLt27d3FfnAFFAVWA1MJMz/9u88MILtG3b1ouIxmQLUdXAC4g8DTwOVANGAI8Ar6nq+NDHyx4REREaHx/vdQxTQKSkQKN6x9g8bxP/YjZTqMueQr/Qa8BumjWLsUEJTZ4hIitVNSLt9IAtDhE5B/gR6AjcBwjwkKpuDElKY/KBmTNh59IdrCOSMJLpzhvUOHc9pUrVx2qGyQ8C9nGoairQR1W/V9UBqvqRFQ1jMnb69Gk++GABdY5PJoxkAMJI5v4Tk0lI8DabMdklmM7xOSLysFj72piA4uPjueWWW/j6695ML1yfJLdBn0RhZhdtRJUq3uYzJrsEUzjaAeOBU3bluDF/d+LECV555RUiIyPZv38/Eyc24/q7biBcltORd4gslkjpyFJERXmd1JjsEbBz3O3juE1Vv8m5SNnPOsdNqCxZsoSYmBg2b95M1apVWb16tTvnzFlVVYAE/M+q6tq1K7GxsTkf1pgsOqvOcVVNFZHeOBcAGmNcR44coXPnzgwYMIBrr72Wr776ivvuu8/rWMbkCOvjMCaLZs+eTXh4OAMHDuTFF19k3bp1VjRMgRLMBYDtgKJAsoicxDklV1X1opAmMyaX+f3332nXrh0jRoygfPnyLFmyhJo1a3ody5gcZ6PjGhOEiRMn8t///pf9+/fTpUsXXnvtNc477zyvYxnjiUwLh4jcmd50VV2U/XGMyV327t1L69atmTBhAlWrVmXWrFlUsfNqTQEXzKGqDn7PzwNqACuBe0OSyJhcQFUZOXIkbdu25fjx47z77ru0b9/eRrE1huAOVT3o/1pEygA9Q5bIGI/t2LGD5s2bM2vWLGrWrElcXBzly5f3OpYxuUYwZ1WltQsIz+4gxngtNTWVgQMHUqFCBRYvXky/fv1YtGiRFQ1j0gimj6M/7pDqOIWmCrAmhJmMyXFbtmyhSZMmLFq0iNq1azN06FDKli3rdSxjcqVgWhzxOH0aK4GlwCuq2jiYjYvI/SKySUS2ikindOaLiPRz568VkWpZWPdlEVERuSyYLMakJzk5mV69elGpUiXWrFlDXFwcc+bMsaJhTAAZtjjcGzaVUNURaaZXEJESqvproA27N38aANTBOby1QkSmquoGv8WigHLuIxIYBERmtq7bz1IH2JGlT2uMn7Vr1xITE0N8fDwNGzZk4MCBXHXVVV7HMibXC9Ti6A+USGd6aaBvENuuAWxV1W2qehoYAzRMs0xD4DN1LAMuEZGSQaz7Ac49QgLfhcqYdJw6dYquXbtSvXp1fvrpJ8aOHcukSZOsaBgTpECFo6KqLkw7UVVnA5WC2HYpYKff613utGCWyXBdEfkX8LOqBuxnEZFmIhIvIvG//hqwcWQKkOXLl1O9enXefPNNHn/8cTZs2MBjjz1md+UzJgsCFY5AJ6wHczJ7ev8npm0hZLRMutNF5AKgC/BGZm+uqkNUNUJVI0qUSK/hZAqS48eP8/LLL3Pbbbdx8OBBpk2bxueff85ll1kXmTFZFahwbBGR+mknikgUsC2Ibe8Cyvi9Lg3sDnKZjKZfD1wLrBGR7e70VSJyZRB5TAG1YMECKlWqRJ8+fWjatCmJiYk88MADXscyJs8KdDruS8A0EXkM54wqgAicIdYbBLHtFUA5EbkW+Bl4AngqzTJTgdYiMganc/yQqu4RkV/TW1dVE4HLz6zsFo8IVd0fRB5TwBw+fJiOHTvy8ccfc/311/P1119zzz33eB3LmDwvw8KhqptFpCLOH/szF/wtBJqr6snMNqyqySLSGpgNFAKGq2qiiLRw5w8GZgD1ga3AceD5QOue5Wc0BdCMGTNo3rw5u3fvpl27drz11ltccMEFXscyJl8IeAfA/MLuAFhw7N+/n7Zt2zJq1CgqVKhAXFwckZGRXscyJk/K6A6AZzPkiDG5jqoybtw4fD4fY8eO5Y033mDlypVWNIwJgWBGxzUmV9uzZw+tWrVi0qRJVK9ena+++opKlYI5Y9wYczaCbnGISJiIVBWRyzNf2pjQU1U++eQTfD4fM2fOpEePHixbtsyKhjEhlmHhEJHBIlLBfX4xzsCGnwGrReTJHMpnTLq2b99OvXr1iI6OpmLFiqxZs4aOHTtSuLA1oo0JtUAtjlp+ZzI9D2xW1YpAdZzhPozJcampqfTv35/w8HCWLl3KgAEDWLBgATfeeKPX0YwpMAL9PDvt97wOMB5AVffa8AzGC5s2bSImJoZvvvmGevXq8fHHH3PNNdd4HcuYAidQi+OgiDQQkarA7cAsABEpDJyfE+GMAUhKSuLdd9+lcuXKbNiwgREjRjBz5kwrGsZ4JFCLoznQD7gSaKuqe93p9wHTQx3MGICEhASio6NZvXo1Dz/8MB999BFXXmkjzBjjpYBXjgP3pzN9Ns4V3caEzMmTJ3n77bfp0aMHl156Kf/73/94+OGHvY5ljCHwjZz6BVpRVV/I/jjGwNKlS4mOjub777/n2Wef5f333+f//u//vI5ljHEFOlTVAlgPjMMZmdZ6xE1IHTt2jC5dutCvXz/KlCnDzJkzuf/+vzV6jTEeC1Q4SgKPAo8DycBYYIKqHsiJYKZgmTdvHk2bNuXHH3+kVatWvPvuu1x44YVexzLGpCPDs6pU9TdVHayq9wDPAZcAiSLynxzKZgqAgwcP0qRJE2rXrk3hwoVZuHAhH330kRUNY3KxTC+zFZFqwJM413LM5M97cxiTJSkpMHMmrF4NVatCUtJUWrduyd69e+nYsSOxsbGcf76d6W1Mbheoc7wbzg2bNgJjgM6qmpxTwUz+EBsbS7du3YBzuIBJlKYUDZlDe+qyCzjOXiCV888/34qGMXlEhvfjEJFUnFvEnnAnnVlQAFXVPDOSnN2Pw3vTpkHXJzex7Gg4YSSTRGGqFlnNW2NuolGjYG5hb4zJaRndjyPQoaprQ5jHFDALFx6m9tFJhOE0WsNI5sGk6SQmhtOokcfhjDFZEqhz/Kf0HsAu4I6ci2jyMlVl6NChDBjQhMnUIcn9rZJEYWYXbUSVKt7mM8ZkXaBh1S8Skc4i8pGI1BVHG5zDV4/lXESTV23bto3atWvTrFkzatT4jTK3XU24LKcj7xBZLJHSkaWIivI6pTEmqwL1cUwBDgBLccanKg4UAV5U1YScCpgdrI8jZ6WkpNC/f3+6dOlCcnIyp0+fGWj5HCAKqAIk4JyklwpA165diY2NzfmwxpgMnU0fx3Xu/TcQkWHAfuBqVT0SoowmH9i4cSPR0dEsW7aM+vXrM3jwYMqUKeN1LGNMNgo0rHrSmSeqmgL8aEXDZCQpKYnu3btTpUoVtmzZwueff860adOsaBiTDwVqcVQWkcPucwHOd1+fOR33opCnM3nCqlWriI6OZs2aNTz22GP079+fyy+3W9Mbk18FOquqkKpe5D4uVNXCfs+taBhOnjxJ586dqVGjBvv27WPSpEmMHTvWioYx+VymQ44Yk54lS5YQExPD5s2biY6Opnfv3hQvXtzrWMaYHBCoj8OYvzl69Cht2rThzjvv5PTp08yZM4e4uDgrGsYUIFY4TNDmzp1LeHg4AwYMoE2bNqxbt446dep4HcsYk8OscJhMHThwgOjoaOrWrct5553H4sWL6du3L8WKFfM6mjHGA1Y4TECTJ0/G5/Px2Wef8eqrr5KQkMDtt9/udSxjjIesc9yk65dffqFNmzaMHz+eKlWqMGPGDKpWrep1LGNMLmAtDvMXqsrnn3+Oz+djypQpdO/eneXLl1vRMMb8wVoc5g87d+6kZcuWTJ8+nZo1axIXF0f58uW9jmWMyWWsxWFITU1l8ODBVKhQgfnz5/Phhx+yaNEiKxrGmHRZi6OA27p1K02aNGHhwoXUrl2bIUOGcO21dg8vY0zGQtriEJH7RWSTiGwVkU7pzBcR6efOXysi1TJbV0R6icj37vKTROSSUH6G/ColJYXevXtTsWJFEhISGDZsGHPmzLGiYYzJVMgKh4gUAgbg3IDBBzwpIr40i0UB5dxHM2BQEOvOBcLde55vBjqH6jPkV+vXr+e2226jQ4cO1K1blw0bNhATE4OIeB3NGJMHhLLFUQPYqqrbVPU0MAZomGaZhsBn6lgGXCIiJQOtq6pzVDXZXX8ZUDqEnyFPS0mBadPgrbecf0+cOE23bt2oVq0a27dvZ8yYMUyePJmrrrrK66jGmDwklH0cpYCdfq93AZFBLFMqyHUBooGx6b25iDTDacVw9dVXZyV3nhYbG0u3bt2Ac7iASZSmFA2ZQ3vqsoufOc6bQCrPPfccjz/+uNdxjTF5UChbHOkd90h7n9qMlsl0XRHpAiQDo9J7c1UdoqoRqhpRokSJIOLmD7GxsagqX36ZQvliN7GeW+nJq6znVkpLGV5/fSmqSs+ePb2OaozJo0JZOHYB/rd/Kw3sDnKZgOuKyLNAA+Bpzeim6QXc6tVQ59hkwnCO6oWRTCPmUKRIDY+TGWPyulAWjhVAORG5VkSKAE8AU9MsMxV4xj276lbgkKruCbSuiNwPvAL8S1WPhzB/nnX48GGWLh3IJK1Nkns0MonCzCnaiCpVvM1mjMn7QlY43A7s1sBsYCMwTlUTRaSFiLRwF5sBbAO2AkOB/wZa113nI+BCYK6IJIjI4FB9hrzo6aef5uKLL2bmzDbs4mcq8B0deYcKfMemo9/z4IOFEBFiY2O9jmqMyaOkIBzpiYiI0Pj4eK9jhNRvv/3GSy+9xMiRI/H5fMTFxXHLLbcycyYkJECVKhAVBYUKeZ3UGJNXiMhKVY1IO92uHM8H/ve//9GqVSt+//13Xn/9dbp06cK5554LQIMGzsMYY7KLFY48bM+ePbRu3ZqJEydSrVo15syZQ+XKlb2OZYzJ52yQwzxIVfn000/x+XxMnz6d9957j++++86KhjEmR1iLI4/56aefaN68ObNnz+aOO+5g2LBh3HTTTV7HMsYUINbiyCNSU1MZMGAA4eHhLFmyhP79+7Nw4UIrGsaYHGctjjxg8+bNxMTEsGTJEurUqcOQIUMoW7as17GMMQWUtThyseTkZHr06EGlSpVYv349n3zyCbNnz7aiYYzxlLU4cqk1a9YQExPDypUradSoEQMGDKBkyZJexzLGGGtx5DanTp3ijTfeICIigp07dzJu3DgmTJhgRcMYk2tYiyMXWbZsGTExMWzYsIHGjRvz4Ycfcumll3odyxhj/sJaHLnAsWPHaNeuHTVr1uTw4cNMnz6dkSNHWtEwxuRK1uLw2Ndff03Tpk3Ztm0bLVq0oEePHlx00UVexzLGmAxZi8Mjhw4dolmzZtx3332ICPPnz2fQoEFWNIwxuZ4VDg9MmzaNChUqEBcXx8svv8zatWu5++67vY5ljDFBscKRg/bv38/TTz/Ngw8+SPHixVm6dCm9evXiggsu8DqaMcYEzQpHDlBVxowZw80338z48ePp2rUrK1eupEYNu42rMSbvsc7xENu9ezctW7Zk6tSp3HLLLcTFxVGxYkWvYxljzFmzFkeIqCpxcXH4fD7mzJlDr169+Pbbb61oGGPyPGtxhMCPP/5I06ZNmTdvHnfeeSfDhg2jXLlyXscyxphsYS2ObJSSkkLfvn0JDw9n+fLlDBw4kPnz51vRMMbkK9biyCYbN24kJiaGpUuXEhUVxccff0yZMmW8jmWMMdnOWhz/UFJSEt27d6dKlSps2rSJkSNHMn36dCsaxph8y1oc/8Dq1auJjo4mISGBRx99lP79+3PFFVd4HcsYY0LKWhxn4eTJk7z66qvccsst7N27lwkTJjBu3DgrGsaYAsFaHFn07bffEhMTw/fff89zzz3H+++/T/Hixb2OZYwxOcYKRwApKTBzJqxeDeXLn2DRos4MGNCPMmXKMGvWLOrVq+d1RGOMyXFWODKQkgKN6h3j5+92UefYZF6nDjv1Xlr+N5X33uvOhRde6HVEY4zxhPVxpBEbG4uIULhwAzbN+55lR8N5TzuxTiMpRWkGDtzGRRddRGxsrNdRjTHGE1Y40oiNjUVVefLJXjRkDmEkAxBGMg/LXN56axqqaoXDGFNgWeFIo0OHDogIo0d3YAp1SXKP5iVRmAlah9dfb4CIWOEwxhRY1sfhUlW++OILPvnkE4oUKUKXLrexfMGNRK5IpN6xScwu2ojykaXYOHsahQp5ndYYY7xjhQPYtWsXLVq0YPr06URGRjJ8+HB8Ph8pXWDmzAtJSHiFN6tAVBRWNIwxBV6BLhypqakMGzaMDh06kJSUxPvvv88LL7xAIbc6FCoEDRo4D2OMMY4CWzh++OEHmjRpwoIFC7j33nsZOnQo1113ndexjDEm1yswhWPaNOdCvkqVUtiypR9vvNGFsLAwhgwZQpMmTRARryMaY0yeIKoauo2L3A/0BQoBw1T1vTTzxZ1fHzgOPKeqqwKtKyL/B4wFygLbgcdU9UCgHJdcVE2v19HUOTaZKVKXHak7ueeBOAYPHkDp0qWz7wMbY0w+IiIrVTUi7fSQnY4rIoWAAUAU4AOeFBFfmsWigHLuoxkwKIh1OwHzVLUcMM99HVDSsdN/XMi3NrUG1597A82bT7aiYYwxZyGU13HUALaq6jZVPQ2MARqmWaYh8Jk6lgGXiEjJTNZtCIxwn48AHsosyIWpB/9yId/9p6bwr3+9johw9913/6MPaYwxBU0o+zhKATv9Xu8CIoNYplQm616hqnsAVHWPiFye3puLSDOcVgzncDHVAQEUSGRkKuzdBhxauHChl/0blwH7vXrzACxX1liurLFcWeNlrmvSmxjKwpHeX+O0HSoZLRPMugGp6hBgCICIxK9M5zid10QkPr3jh16zXFljubLGcmVNbswVykNVuwD/+6eWBnYHuUygdX9xD2fh/rsvGzMbY4zJRCgLxwqgnIhcKyJFgCeAqWmWmQo8I45bgUPuYahA604FnnWfPwtMCeFnMMYYk0bIDlWparKItAZm45xSO1xVE0WkhTt/MDAD51TcrTin4z4faF130+8B40QkBtgBPBpEnCHZ98myleXKGsuVNZYrayxXkEJ6HYcxxpj8x4ZVN8YYkyVWOIwxxmRJni4cInK/iGwSka0i8rcryN1O937u/LUiUi3YdUOc62k3z1oR+VZEKvvN2y4i60QkQUTiczjX3SJyyH3vBBF5I9h1cyBbB79c60UkxR1+JmTfmYgMF5F9IrI+g/le7V+Z5fJq/8oslyf7VxC5cnzfcrddRkTmi8hGEUkUkRfTWcaTfSxTqponHzid5j8A1wFFgDWAL80y9YGZONeF3Ap8F+y6Ic5VEyjuPo86k8t9vR24zKPv625g2tmsG+psaZZ/EPg6B76zO4FqwPoM5uf4/hVkrhzfv4LM5dX+FTCXF/uWu+2SQDX3+YXA5tzwNyyYR15ucYRqSJOQ51LVb/XPgRmX4VynEmr/5DOH8vs6m+0/CYzOxvdPl6ouAn4PsIgX+1emuTzav4L5vjLi6feVRo7sW+CMfKHuoK6qegTYiDNqhj9P9rHM5OXCkdFwJcEsE8y6oczlLwbnF8UZCswRkZXiDJuSXYLNdZuIrBGRmSJSIYvrhjobInIBcD8wwW9yqL6zzHixf2VVTu1fwfJi/wqKl/uWiJQFqgLfpZmVK/exvHw/Dk+HNAkg6G2LyD04/2Pf4Tf5dlXdLc4YXHNF5Hv3F1NO5FoFXKOqR0WkPjAZZ+TiUH5fwWY740HgG1X1/wUZqu8sM17sX0HL4f0rGF7tX8HyZN8SkWI4xaqtqh5OOzudVTzfx/JyiyNUQ5rkRC5EpBIwDGioqr+dma6qu91/9wGTcJqkOZJLVQ+r6lH3+QwgTEQuC2bdUGfz8wRpDiWE8DvLjBf7V1A82L8y5eH+Fawc37dEJAynaIxS1YnpLJI797Gc6kzJ7gdOa2kbcC1/dg5VSLPMA/y1Y2l5sOuGONfVOFfL10wzvShwod/zb4H7czDXlfx5UWgNnCvzJZTfV1b+ewAX4xyrLpoT35m7zbJk3Nmb4/tXkLlyfP8KMpcn+1dmuTzctwT4DPgwwDKe7WOBHnn2UJWGbkiTnMj1BnApMFCcId2T1Rn98gpgkjutMPCFqs7KwVyPAC1FJBk4ATyhzl4asu8rC9kAGgFzVPWY3+oh+85EZDTOmUCXicguoCsQ5pcpx/evIHPl+P4VZC5P9q8gckEO71uu24H/AOtEJMGd9ipO4fd0H8uMDTlijDEmS/JyH4cxxhgPWOEwxhiTJVY4jDHGZIkVDmOMMVlihcMYY0yWWOEweZ47mumZkU2/FJFL3OllRURF5C2/ZS8TkSQR+ch9fZOILHDX3ygif7vbmoj8KCI3pZn2oYh0DJDpqF+Gp7LpowYkIm1F5JmzXHeBiES4z18NYvkxIlLubN7L5H1WOEx+cEJVq6hqOM5FXK385m0DGvi9fhTwP9+9H/CBu/7NQP90tj8G56piAETkHJxrEsYGka0sEPLCISKFgWjgi2zYXKaFAxgEZFg4Tf5mhcPkN0v562BvJ4CNZ35NA48D4/zml8QZvgEAVV2XzjZH41c4cIbp3q6qP4lIO7els15E2qaz7ntALbdF85LbAlksIqvcR01wipGIDBTnvgzTRGSGiDzizqsuIgvdgfZmu6OjpnUvsMq9MOxmEVl+Zob7nmvd5/eJyGpx7jExXETO9d+IiLwHnO/mHSUiRUVkujsw4XoRedxddDFQ2y1YpoCxwmHyDREpBNwHTE0zawzwhIiUBlL465g+HwBfizNa60tnDnP5U9W1QKr8eUOkJ4DRIlId50reSJzhIJqKSNU0q3cCFrstmg+AfUAdVa2GU8T6ucv9G6d1UhFoAtzmfqYwnFbQI6paHRgOdE/n498OrHTzbgSKiMh17rzHgXEich7wKfC4qlbEuRq6ZZrP2ok/W3BP44wWu1tVK7stulnucqk4VzNXxhQ4VjhMfnC+O2TDb8D/AXPTzJ8F1MG518JfDi+p6ifAzcB4nGEplqX9Fe4ajVN8CuPc92A8zqizk1T1mDqD900EamWSNQwYKiLr3G343Ol3AONVNVVV9wLz3ek3AeE4I7MmAK+R/v01SgK/+r0eBzzmPn/c/dw3AT+q6mZ3+gic1lMg63BaFj1EpJaqHvKbtw+4KpP1TT5khcPkBydUtQpwDc6Ab/59HKhzo5uVQHv+eq+FM/N3q+pwVW0IJOP8oU5rNM4f4trAWnVGS01vaOvMvAT8gvNLPcLNS4BtCZDotgCqqGpFVa2bznIngPP8Xo8FHhORGwFV1S1nk9ctMtVxCsi74ne7V/f9TmR1mybvs8Jh8g331/ALwMvuIR5/fYBX1G+Icfjjvs1h7vMrcQYH/Dmdbf+A06J5jz+H3l4EPCQiF4hIUZyB8hanWfUIzm1Bz7gY2OMe6vkPzgB1AEuAh92+jitwWj8Am4ASIvLHoSv58wZI/jYCN6TJmwK8zp+trO+BsiJyZrn/AAvT2VaS33dyFXBcVT8HeuPcgvWMG/nriQamgLCOLZOvqOpqEVmD0w+x2G96Iun/kasL9BWRk+7rDu6hovSMBt7FuS8DqrpKRD4FznRED1PV1WnWWYsz+usanP6FgcAEEXkU53DUmdFYJ+D0z6zHuff0d8AhVT3tdpL3E5GLcf6f/TCdzzITGJlm2ligF87Q26jqSRF5HhjvHnJbAQzm74YAa0VkFc6w371EJBVIwu0TcYvbCVXdk8F3ZfIxGx3XmFxCRIqpc3e8S3GK0e0Bilh6608COrqHpUJKRF4CDqtqXKjfy+Q+1uIwJveY5p7VVQR4KytFw9UJp5M85IUDOMjfWzimgLAWhzHGmCyxznFjjDFZYoXDGGNMlljhMMYYkyVWOIwxxmSJFQ5jjDFZ8v8icOSSAWOHpwAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plt.errorbar(Vrms, Irms, errIrms, fmt = 'ko', markersize = 5,\\\n", " linewidth = 1.8,\\\n", " markeredgecolor = 'b',\\\n", " markerfacecolor = 'r',\\\n", " capsize = 5)\n", "\n", "# Plot the best-fit line...\n", "xx = np.arange(-1, 10, 0.1)\n", "plt.plot(xx, fitFcn(xx), 'k-')\n", "\n", "# Do some formatting...\n", "plt.xlabel('RMS Voltage (volts)')\n", "plt.ylabel('RMS Current (amps)')\n", "plt.title('Weighted Linear Fit')\n", "plt.axis((0, 2.2, 0, 0.01));" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Save the plot as a pdf." ] }, { "cell_type": "code", "execution_count": 89, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.savefig('Weighted linear fit plot.pdf', format='pdf')" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.8.13" } }, "nbformat": 4, "nbformat_minor": 4 }